26 #ifndef O2SCL_PROB_DENS_FUNC_H
27 #define O2SCL_PROB_DENS_FUNC_H
29 #include <gsl/gsl_rng.h>
30 #include <gsl/gsl_randist.h>
31 #include <gsl/gsl_cdf.h>
33 #include <boost/numeric/ublas/vector.hpp>
35 #include <o2scl/hist.h>
36 #include <o2scl/rng_gsl.h>
37 #include <o2scl/search_vec.h>
39 #ifndef DOXYGEN_NO_O2NS
55 virtual double sample()
const=0;
61 virtual double cdf(
double x)
const=0;
67 virtual double entropy()
const=0;
105 r=gsl_rng_alloc(gsl_rng_mt19937);
114 O2SCL_ERR2(
"Tried to create a Gaussian dist. with sigma",
115 "<0 in prob_dens_gaussian::prob_dens_gaussian().",
120 r=gsl_rng_alloc(gsl_rng_mt19937);
131 r=gsl_rng_alloc(gsl_rng_mt19937);
137 if (
this==&pdg)
return *
this;
157 "in prob_dens_gaussian::prob_dens_gaussian().",
171 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
180 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
189 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
196 virtual double cdf(
double x)
const {
198 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
207 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
210 if (in_cdf<0.0 || in_cdf>1.0) {
211 O2SCL_ERR2(
"Requested cdf inverse outside of [0,1] in ",
212 "prob_dens_gaussian::invert_cdf().",
exc_einval);
220 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
274 r=gsl_rng_alloc(gsl_rng_mt19937);
288 r=gsl_rng_alloc(gsl_rng_mt19937);
299 r=gsl_rng_alloc(gsl_rng_mt19937);
305 if (
this==&pdg)
return *
this;
333 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
342 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
351 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
354 return gsl_ran_flat(
r,
ll,
ul);
360 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
363 if (x<ll || x>
ul)
return 0.0;
364 return gsl_ran_flat_pdf(x,
ll,ul);
368 virtual double cdf(
double x)
const {
370 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
373 if (x<
ll)
return 0.0;
374 if (x>
ul)
return 1.0;
375 return gsl_cdf_flat_P(x,
ll,
ul);
381 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
384 if (in_cdf<0.0 || in_cdf>1.0) {
385 O2SCL_ERR2(
"Requested cdf inverse outside of [0,1] in ",
386 "prob_dens_uniform::invert_cdf().",
exc_einval);
388 return gsl_cdf_flat_Pinv(in_cdf,
ll,
ul);
394 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
447 r=gsl_rng_alloc(gsl_rng_mt19937);
457 O2SCL_ERR2(
"Tried to create log normal dist. with mu or sigma",
458 "<0 in prob_dens_lognormal::prob_dens_lognormal().",
463 r=gsl_rng_alloc(gsl_rng_mt19937);
474 r=gsl_rng_alloc(gsl_rng_mt19937);
480 if (
this==&pdg)
return *
this;
490 O2SCL_ERR2(
"Tried to set mu or sigma negative",
491 "in prob_dens_lognormal::prob_dens_lognormal().",
513 return gsl_ran_lognormal_pdf(x,
mu_,
sigma_);
517 virtual double cdf(
double x)
const {
521 return gsl_cdf_lognormal_P(x,
mu_,
sigma_);
526 if (in_cdf<0.0 || in_cdf>1.0) {
527 O2SCL_ERR2(
"Requested cdf inverse outside of [0,1] in ",
528 "prob_dens_lognormal::invert_cdf().",
exc_einval);
530 return gsl_cdf_lognormal_Pinv(in_cdf,
mu_,
sigma_);
536 O2SCL_ERR2(
"Parameters not set in prob_dens_lognormal::",
587 virtual double sample()
const;
599 virtual double cdf(
double x)
const;
611 #ifndef DOXYGEN_NO_O2NS
double sigma_
Width parameter.
virtual double entropy() const
The inverse cumulative distribution function.
virtual double sample() const
Sample from the specified density.
prob_dens_gaussian & operator=(const prob_dens_gaussian &pdg)
Copy constructor with operator=.
virtual double lower_limit() const
Lower limit of the range.
void set_sigma(double sigma)
Set the Gaussian width.
invalid argument supplied by user
ubvector sum
Normalized partial sum of histogram bins.
size_t n
Number of original histogram bins.
prob_dens_lognormal & operator=(const prob_dens_lognormal &pdg)
Copy constructor with operator=.
prob_dens_lognormal()
Create a blank lognormal distribution.
A one-dimensional Gaussian probability density.
prob_dens_gaussian()
Create a standard normal distribution.
void set_mu_sigma(double mu, double sigma)
Set the maximum and width of the lognormal distribution.
A one-dimensional histogram class.
prob_dens_gaussian(const prob_dens_gaussian &pdg)
Copy constructor.
prob_dens_lognormal(const prob_dens_lognormal &pdg)
Copy constructor.
double get_sigma()
Get the Gaussian width.
virtual double sample() const
Sample from the specified density.
prob_dens_gaussian(double cent, double sigma)
Create a Gaussian distribution with width sigma.
virtual double entropy() const
The inverse cumulative distribution function.
virtual double operator()(double x) const
The normalized density.
virtual double cdf(double x) const
The cumulative distribution function (from the lower tail)
virtual double entropy() const =0
Entropy of the distribution ( )
virtual double cdf(double x) const
The cumulative distribution function (from the lower tail)
A one-dimensional probability density function.
A one-dimensional probability density over a finite range.
virtual double upper_limit() const =0
Uower limit of the range.
double get_center()
Get the center.
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
virtual double operator()(double x) const =0
The normalized density.
virtual double operator()(double x) const
The normalized density.
virtual double invert_cdf(double x) const
Inverse cumulative distribution function (from the lower tail)
ubvector range
Vector specifying original histogram bins.
double cent_
Central value.
search_vec< ubvector > sv
Search through the partial sums.
virtual double invert_cdf(double in_cdf) const
The inverse cumulative distribution function.
double sigma_
Width parameter.
virtual double upper_limit() const
Uower limit of the range.
A one-dimensional probability density over the positive real numbers.
virtual double operator()(double x) const
The normalized density.
virtual double cdf(double x) const =0
The cumulative distribution function (from the lower tail)
Random number generator (GSL)
Searching class for monotonic data with caching.
prob_dens_lognormal(double mu, double sigma)
Create lognormal distribution with mean parameter mu and width parameter sigma.
virtual double sample() const
Generate a sample.
virtual double invert_cdf(double in_cdf) const
The inverse cumulative distribution function.
rng_gsl rng
Random number generator.
gsl_rng * r
Base GSL random number generator.
Probability density function based on a histogram.
void init(hist &h)
Initialize with histogram h.
virtual double lower_limit() const =0
Lower limit of the range.
void set_seed(unsigned long int s)
Set the seed.
virtual double sample() const =0
Sample from the specified density.
void set_seed(unsigned long int s)
Set the seed.
virtual double cdf(double x) const
Cumulative distribution function (from the lower tail)
virtual double entropy() const
Inverse cumulative distribution function (from the lower tail)
virtual double invert_cdf(double cdf) const =0
The inverse cumulative distribution function.
Lognormal density function.
void set_center(double cent)
Set the center.
gsl_rng * r
The GSL random number generator.