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> 34 #include <gsl/gsl_poly.h> 38 #include <boost/numeric/ublas/vector.hpp> 39 #include <boost/numeric/ublas/operation.hpp> 41 #include <o2scl/hist.h> 42 #include <o2scl/rng_gsl.h> 43 #include <o2scl/search_vec.h> 44 #include <o2scl/cholesky.h> 47 #ifndef DOXYGEN_NO_O2NS 74 virtual double pdf(
double x)
const {
84 virtual double cdf(
double x)
const {
142 O2SCL_ERR2(
"Tried to create a Gaussian dist. with sigma",
143 "<0 in prob_dens_gaussian::prob_dens_gaussian().",
187 "in prob_dens_gaussian::prob_dens_gaussian().",
202 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
211 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
214 return cent_+gsl_ran_gaussian(&r,sigma_);
218 virtual double pdf(
double x)
const {
220 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
223 return gsl_ran_gaussian_pdf(x-cent_,sigma_);
229 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
232 return log(gsl_ran_gaussian_pdf(x-cent_,sigma_));
236 virtual double cdf(
double x)
const {
238 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
241 return gsl_cdf_gaussian_P(x-cent_,sigma_);
247 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
250 if (in_cdf<0.0 || in_cdf>1.0) {
251 O2SCL_ERR2(
"Requested cdf inverse outside of [0,1] in ",
252 "prob_dens_gaussian::invert_cdf().",
exc_einval);
254 return gsl_cdf_gaussian_Pinv(in_cdf,sigma_)+cent_;
260 O2SCL_ERR2(
"Width not set in prob_dens_gaussian::",
278 virtual double lower_limit()
const=0;
281 virtual double upper_limit()
const=0;
341 if (
this==&pdg)
return *
this;
370 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
379 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
388 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
391 return gsl_ran_flat(&r,ll,ul);
395 virtual double pdf(
double x)
const {
397 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
400 if (x<ll || x>ul)
return 0.0;
401 return gsl_ran_flat_pdf(x,ll,ul);
407 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
410 if (x<ll || x>ul)
return 0.0;
411 return log(gsl_ran_flat_pdf(x,ll,ul));
415 virtual double cdf(
double x)
const {
417 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
420 if (x<ll)
return 0.0;
421 if (x>ul)
return 1.0;
422 return gsl_cdf_flat_P(x,ll,ul);
428 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
431 if (in_cdf<0.0 || in_cdf>1.0) {
432 O2SCL_ERR2(
"Requested cdf inverse outside of [0,1] in ",
433 "prob_dens_uniform::invert_cdf().",
exc_einval);
435 return gsl_cdf_flat_Pinv(in_cdf,ll,ul);
441 O2SCL_ERR2(
"Limits not set in prob_dens_uniform::",
503 O2SCL_ERR2(
"Tried to create log normal dist. with mu or sigma",
504 "<0 in prob_dens_lognormal::prob_dens_lognormal().",
523 if (
this==&pdg)
return *
this;
533 O2SCL_ERR2(
"Tried to set mu or sigma negative",
534 "in prob_dens_lognormal::prob_dens_lognormal().",
548 return gsl_ran_lognormal(&r,mu_,sigma_);
552 virtual double pdf(
double x)
const {
556 return gsl_ran_lognormal_pdf(x,mu_,sigma_);
564 return log(gsl_ran_lognormal_pdf(x,mu_,sigma_));
568 virtual double cdf(
double x)
const {
572 return gsl_cdf_lognormal_P(x,mu_,sigma_);
577 if (in_cdf<0.0 || in_cdf>1.0) {
578 O2SCL_ERR2(
"Requested cdf inverse outside of [0,1] in ",
579 "prob_dens_lognormal::invert_cdf().",
exc_einval);
581 return gsl_cdf_lognormal_Pinv(in_cdf,mu_,sigma_);
587 O2SCL_ERR2(
"Parameters not set in prob_dens_lognormal::",
647 virtual double lower_limit()
const;
650 virtual double upper_limit()
const;
653 virtual double pdf(
double x)
const;
656 virtual double log_pdf(
double x)
const;
659 virtual double cdf(
double x)
const;
675 template<
class vec_t=boost::numeric::ublas::vector<
double> >
681 virtual size_t dim()
const {
686 virtual double pdf(
const vec_t &x)
const {
691 virtual double log_pdf(
const vec_t &x)
const {
705 template<
class vec_t=boost::numeric::ublas::vector<
double> >
711 std::vector<prob_dens_func>
list;
720 virtual size_t dim()
const {
725 virtual double pdf(
const vec_t &x)
const {
727 for(
size_t i=0;i<list.size();i++) ret*=list[i].
pdf(x[i]);
732 virtual double log_pdf(
const vec_t &x)
const {
734 for(
size_t i=0;i<list.size();i++) ret*=list[i].
pdf(x[i]);
740 for(
size_t i=0;i<list.size();i++) x[i]=list[i]();
798 template<
class vec_t=boost::numeric::ublas::vector<
double> >
828 void set(
double x_cent,
double y_cent,
double x_std,
double y_std,
830 if (fabs(covar)>=1.0) {
831 O2SCL_ERR2(
"Covariance cannot have magnitude equal or larger than ",
832 "1 in prob_dens_mdim_biv_gaussian::set().",
845 virtual double pdf(
const vec_t &v)
const {
846 double x=v[0], y=v[1];
847 double arg=-((x-x0)*(x-x0)/sig_x/sig_x+
848 (y-y0)*(y-y0)/sig_y/sig_y-
849 2.0*rho*(x-x0)*(y-y0)/sig_x/sig_y)/2.0/(1.0-rho*rho);
862 double arg=-2.0*log(1.0-integral);
866 sig_y/sqrt(1.0-rho*rho);
872 virtual void contour(
double level,
double theta, vec_t &x) {
874 O2SCL_ERR2(
"Cannot produce contours for negative values in ",
875 "prob_dens_mdim_biv_gaussian::contour().",
880 O2SCL_ERR2(
"Cannot produce contours larger than maximum in ",
881 "prob_dens_mdim_biv_gaussian::contour().",
885 sqrt(1.0-rho*rho))*2.0*(1.0-rho*rho);
886 double r2=arg/(cos(theta)*cos(theta)/sig_x/sig_x+
887 sin(theta)*sin(theta)/sig_y/sig_y-
888 2.0*rho/sig_x/sig_y*cos(theta)*sin(theta));
889 x[0]=sqrt(r2)*cos(theta)+x0;
890 x[1]=sqrt(r2)*sin(theta)+y0;
922 template<
class vec_t=boost::numeric::ublas::vector<
double>,
923 class mat_t=boost::numeric::ublas::matrix<
double> >
955 virtual size_t dim()
const {
967 set(p_ndim,p_peak,covar);
972 void set(
size_t p_ndim, vec_t &p_peak, mat_t &covar) {
974 O2SCL_ERR(
"Zero dimension in prob_dens_mdim_gaussian::set().",
980 for(
size_t i=0;i<ndim;i++) peak[i]=p_peak[i];
990 o2scl_linalg::cholesky_invert<mat_t>(ndim,covar_inv);
994 for(
size_t i=0;i<ndim;i++) {
996 for(
size_t j=0;j<ndim;j++) {
997 if (i<j) chol(i,j)=0.0;
1009 void set_alt(
size_t p_ndim, vec_t &p_peak, mat_t &p_chol,
1010 mat_t &p_covar_inv,
double p_norm) {
1014 covar_inv=p_covar_inv;
1023 template<
class vec_vec_t,
class mat_col_t,
class func_t>
1025 vec_vec_t &x, vec_t &y, func_t &fcovar) {
1029 mat_t KXsX(n_dim,n_init);
1030 for(
size_t irow=n_init;irow<n_dim+n_init;irow++) {
1031 for(
size_t icol=0;icol<n_init;icol++) {
1032 KXsX(irow-n_init,icol)=fcovar(x[irow],x[icol]);
1036 mat_t KXXs=boost::numeric::ublas::trans(KXsX);
1038 mat_t KXX(n_init,n_init);
1039 for(
size_t irow=0;irow<n_init;irow++) {
1040 for(
size_t icol=0;icol<n_init;icol++) {
1042 KXX(irow,icol)=KXX(icol,irow);
1044 KXX(irow,icol)=fcovar(x[irow],x[icol]);
1049 mat_t KXsXs(n_dim,n_dim);
1050 for(
size_t irow=n_init;irow<n_dim+n_init;irow++) {
1051 for(
size_t icol=n_init;icol<n_dim+n_init;icol++) {
1053 KXsXs(irow-n_init,icol-n_init)=KXsXs(icol-n_init,irow-n_init);
1055 KXsXs(irow-n_init,icol-n_init)=fcovar(x[irow],x[icol]);
1061 mat_t inv_KXX(n_init,n_init);
1067 "prob_dens_mdim_gaussian::set_gproc().",
1070 o2scl_linalg::LU_invert<mat_t,mat_t,mat_col_t>(n_init,KXX,p,inv_KXX);
1073 vec_t prod(n_init), mean(n_dim);
1074 boost::numeric::ublas::axpy_prod(inv_KXX,y,prod,
true);
1075 boost::numeric::ublas::axpy_prod(KXsX,prod,mean,
true);
1078 mat_t covar(n_dim,n_dim), prod2(n_init,n_dim), prod3(n_dim,n_dim);
1079 boost::numeric::ublas::axpy_prod(inv_KXX,KXXs,prod2,
true);
1080 boost::numeric::ublas::axpy_prod(KXsX,prod2,prod3,
true);
1084 this->
set(n_dim,mean,covar);
1089 virtual double pdf(
const vec_t &x)
const {
1091 O2SCL_ERR2(
"Distribution not set in prob_dens_mdim_gaussian::",
1095 for(
size_t i=0;i<ndim;i++) q[i]=x[i]-peak[i];
1096 vtmp=prod(covar_inv,q);
1097 ret*=exp(-0.5*inner_prod(q,vtmp));
1104 O2SCL_ERR2(
"Distribution not set in prob_dens_mdim_gaussian::",
1107 double ret=log(norm);
1108 for(
size_t i=0;i<ndim;i++) q[i]=x[i]-peak[i];
1109 vtmp=prod(covar_inv,q);
1110 ret+=-0.5*inner_prod(q,vtmp);
1117 O2SCL_ERR2(
"Distribution not set in prob_dens_mdim_gaussian::",
1120 for(
size_t i=0;i<ndim;i++) q[i]=pdg();
1122 for(
size_t i=0;i<ndim;i++) x[i]=peak[i]+vtmp[i];
1145 template<
class vec_t=boost::numeric::ublas::vector<
double> >
1156 virtual double pdf(
const vec_t &x_B,
const vec_t &x_A)
const=0;
1159 virtual double log_pdf(
const vec_t &x_B,
const vec_t &x_A)
const=0;
1162 virtual void operator()(
const vec_t &x_B, vec_t &x_A)
const=0;
1236 template<
class vec_t=boost::numeric::ublas::vector<
double> >
1268 for(
size_t i=0;i<sz;i++) {
1271 if (!std::isfinite(low[i]) || !std::isfinite(high[i])) {
1272 O2SCL_ERR2(
"Limit not finite in prob_cond_mdim_fixed_step::",
1277 if (low[i]>high[i]) {
1304 (vec_t &step, vec_t &low, vec_t &high) {
1305 if (step.size()!=low.size()) {
1306 O2SCL_ERR2(
"Vectors 'step' and 'low' mismatched in ",
1307 "prob_cond_mdim_fixed_step constructor.",
1310 if (step.size()!=high.size()) {
1311 O2SCL_ERR2(
"Vectors 'step' and 'high' mismatched in ",
1312 "prob_cond_mdim_fixed_step constructor.",
1315 set_internal(step.size(),step,low,high);
1320 virtual int set(vec_t &step, vec_t &low, vec_t &high) {
1321 if (step.size()!=low.size()) {
1322 O2SCL_ERR2(
"Vectors 'step' and 'low' mismatched in ",
1323 "prob_cond_mdim_fixed_step::set().",
1326 if (step.size()!=high.size()) {
1327 O2SCL_ERR2(
"Vectors 'step' and 'high' mismatched in ",
1328 "prob_cond_mdim_fixed_step::set().",
1331 set_internal(step.size(),step,low,high);
1337 return u_step.size();
1341 virtual double pdf(
const vec_t &x,
const vec_t &
x2)
const {
1343 for(
size_t i=0;i<u_step.size();i++) {
1345 if (fabs(x2[i]-x[i])>u_step[i])
return 0.0;
1347 if (x[i]-u_step[i]<u_low[i]) {
1348 if (x[i]+u_step[i]>u_high[i]) {
1350 vol1*=u_high[i]-u_low[i];
1353 vol1*=x[i]+u_step[i]-u_low[i];
1356 if (x[i]+u_step[i]>u_high[i]) {
1358 vol1*=u_high[i]-(x[i]-u_step[i]);
1362 vol1*=2.0*u_step[i];
1370 virtual double log_pdf(
const vec_t &x,
const vec_t &
x2)
const {
1371 return log(
pdf(x,x2));
1376 size_t nv=u_step.size();
1377 for(
size_t i=0;i<nv;i++) {
1378 if (x[i]<u_low[i] || x[i]>u_high[i]) {
1379 std::cout <<
"Herex: " << i <<
" " << x[i] <<
" " 1380 << u_low[i] <<
" " << u_high[i] << std::endl;
1381 O2SCL_ERR(
"Input out of bounds in fixed_step::operator().",
1385 x2[i]=x[i]+u_step[i]*(rg.
random()*2.0-1.0);
1386 }
while (x2[i]<u_low[i] || x2[i]>u_high[i]);
1404 template<
class vec_t=boost::numeric::ublas::vector<
double> >
1425 virtual double pdf(
const vec_t &x,
const vec_t &
x2)
const {
1426 return base.
pdf(x2);
1430 virtual double log_pdf(
const vec_t &x,
const vec_t &
x2)
const {
1449 template<
class vec_t=boost::numeric::ublas::vector<
double>,
1450 class mat_t=boost::numeric::ublas::matrix<
double> >
1497 void set(
size_t p_ndim, mat_t &covar) {
1499 O2SCL_ERR(
"Zero dimension in prob_cond_mdim_gaussian::set().",
1513 o2scl_linalg::cholesky_invert<mat_t>(ndim,covar_inv);
1517 for(
size_t i=0;i<ndim;i++) {
1519 for(
size_t j=0;j<ndim;j++) {
1520 if (i<j) chol(i,j)=0.0;
1530 virtual double pdf(
const vec_t &x,
const vec_t &
x2)
const {
1532 O2SCL_ERR2(
"Distribution not set in prob_cond_mdim_gaussian::",
1536 for(
size_t i=0;i<ndim;i++) q[i]=x2[i]-x[i];
1537 vtmp=prod(covar_inv,q);
1538 ret*=exp(-0.5*inner_prod(q,vtmp));
1543 virtual double log_pdf(
const vec_t &x,
const vec_t &
x2)
const {
1545 O2SCL_ERR2(
"Distribution not set in prob_cond_mdim_gaussian::",
1548 double ret=log(norm);
1549 for(
size_t i=0;i<ndim;i++) q[i]=x2[i]-x[i];
1550 vtmp=prod(covar_inv,q);
1551 ret+=-0.5*inner_prod(q,vtmp);
1558 O2SCL_ERR2(
"Distribution not set in prob_cond_mdim_gaussian::",
1561 for(
size_t i=0;i<ndim;i++) q[i]=pdg();
1563 for(
size_t i=0;i<ndim;i++) x2[i]=x[i]+vtmp[i];
1569 #ifdef O2SCL_NEVER_DEFINED 1576 template<
class vec_t=boost::numeric::ublas::vector<
double>,
1577 class mat_t=boost::numeric::ublas::matrix<
double>,
1578 class mat_col_t=boost::numeric::ublas::matrix_column<mat_t> >
1579 class prob_dens_mdim_gproc :
1587 #ifndef DOXYGEN_NO_O2NS vec_t peak
Location of the peak.
double sigma_
Width parameter.
virtual double entropy() const
The inverse cumulative distribution function.
void set_gproc(size_t n_dim, size_t n_init, vec_vec_t &x, vec_t &y, func_t &fcovar)
Given a data set and a covariance function, construct probability distribution based on a Gaussian pr...
rng_gsl rg
Internal random number generator.
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
virtual double cdf(double x) const
The cumulative distribution function (from the lower tail)
prob_cond_mdim_gaussian(size_t p_ndim, mat_t &covar)
Create a distribution from the covariance matrix.
double mean()
Get the center.
virtual double log_pdf(const vec_t &x) const
The log of the normalized density.
virtual double invert_cdf(double in_cdf) const
The inverse cumulative distribution function.
virtual size_t dim() const
The dimensionality.
int set_internal(size_t sz, vec_t &step, vec_t &low, vec_t &high)
Internal set function.
virtual void operator()(const vec_t &x, vec_t &x2) const
Sample the distribution.
virtual void operator()(const vec_t &x, vec_t &x2) const
Sample the distribution.
virtual double pdf(const vec_t &v) const
Compute the normalized probability density.
virtual double pdf(const vec_t &x) const
The normalized density.
virtual size_t dim() const
The dimensionality.
Conditional probability for a random walk inside a hypercube.
void set_seed(unsigned long int s)
Set the seed.
prob_dens_gaussian & operator=(const prob_dens_gaussian &pdg)
Copy constructor with operator=.
o2scl::prob_dens_gaussian pdg
Standard normal.
virtual double log_pdf(const vec_t &x, const vec_t &x2) const
The log of the conditional probability.
void set_sigma(double sigma)
Set the Gaussian width (must be positive)
invalid argument supplied by user
prob_dens_mdim_gaussian(size_t p_ndim, vec_t &p_peak, mat_t &covar)
Create a distribution from the covariance matrix.
ubvector sum
Normalized partial sum of histogram bins.
A multi-dimensional conditional probability density function independent of the input.
size_t n
Number of original histogram bins.
vec_t vtmp
Temporary storage 2.
virtual double entropy() const
Inverse cumulative distribution function (from the lower tail)
const ubvector & bin_ranges()
Get reference to bin ranges.
virtual double level_fixed_integral(double integral)
Return the contour level corresponding to a fixed integral.
prob_dens_lognormal & operator=(const prob_dens_lognormal &pdg)
Copy constructor with operator=.
virtual void operator()(vec_t &x) const
Sample the distribution.
A class for representing permutations.
double sig_x
The x standard deviation.
prob_dens_lognormal()
Create a blank lognormal distribution.
int LU_decomp(const size_t N, mat_t &A, o2scl::permutation &p, int &signum)
Compute the LU decomposition of the matrix A.
o2scl::rng_gsl r
Base GSL random number generator.
o2scl::prob_dens_gaussian pdg
Standard normal.
double random()
Return a random number in .
double x0
The x coordinate of the centroid.
virtual double pdf(const vec_t &x, const vec_t &x2) const
The conditional probability.
A one-dimensional Gaussian probability density.
A multi-dimensional conditional probability density function.
virtual double cdf(double x) const
The cumulative distribution function (from the lower tail)
rng_gsl r
The GSL random number generator.
prob_dens_gaussian()
Create a standard normal distribution.
mat_t chol
Cholesky decomposition.
virtual size_t dim() const
The dimensionality.
virtual double log_pdf(const vec_t &x) const
The log of the normalized density.
mat_t covar_inv
Inverse of the covariance matrix.
virtual double invert_cdf(double cdf) const
The inverse cumulative distribution function.
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.
virtual double pdf(const vec_t &x) const
The normalized density.
prob_dens_lognormal(const prob_dens_lognormal &pdg)
Copy constructor.
std::vector< double > u_step
Step sizes.
vec_t vtmp
Temporary storage 2.
virtual double invert_cdf(double in_cdf) const
The inverse cumulative distribution function.
virtual size_t dim() const
The dimensionality.
prob_cond_mdim_gaussian()
Create an empty distribution.
prob_dens_gaussian(double cent, double sigma)
Create a Gaussian distribution with width sigma.
A multi-dimensional Gaussian conditional probability density function.
virtual double log_pdf(const vec_t &x, const vec_t &x2) const
The log of the conditional probability.
virtual double log_metrop_hast(const vec_t &x_B, vec_t &x_A) const
Sample the distribution and return the log of the Metropolis-Hastings ratio.
A multidimensional distribution formed by the product of several one-dimensional distributions.
double norm
Normalization factor.
A one-dimensional probability density function.
virtual double pdf(double x) const
The normalized density.
virtual double entropy() const
The inverse cumulative distribution function.
A one-dimensional probability density over a finite range.
virtual double entropy() const
Entropy of the distribution ( )
virtual void operator()(vec_t &x) const
Sample the distribution.
void cholesky_decomp(const size_t M, mat_t &A)
Compute the in-place Cholesky decomposition of a symmetric positive-definite square matrix...
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
A multi-dimensional probability density function.
void set_seed(unsigned long int s)
Set the random number generator seed.
virtual double cdf(double x) const
The cumulative distribution function (from the lower tail)
ubvector range
Vector specifying original histogram bins.
A bivariate gaussian probability distribution.
double cent_
Central value.
virtual double pdf(const vec_t &x) const
The normalized density.
virtual size_t dim() const
Return the dimensionality.
search_vec< ubvector > sv
Search through the partial sums.
virtual double log_pdf(double x) const
The log of the normalized density.
virtual double operator()() const
Sample from the specified density.
std::vector< double > u_high
Upper limits.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
virtual double operator()() const
Sample from the specified density.
double sigma_
Width parameter.
virtual double pdf(const vec_t &x, const vec_t &x2) const
The conditional probability.
vec_t q
Temporary storage 1.
virtual double log_pdf(double x) const
The log of the normalized density.
double rho
The covariance.
A one-dimensional probability density over the positive real numbers.
vec_t q
Temporary storage 1.
virtual void operator()(const vec_t &x, vec_t &x2) const
Sample the distribution.
virtual double log_pdf(const vec_t &x) const
The log of the normalized density.
A multi-dimensional Gaussian probability density function using a Cholesky decomposition.
virtual double pdf(const vec_t &x, const vec_t &x2) const
The conditional probability.
virtual size_t dim() const
The dimensionality.
prob_cond_mdim_indep(prob_dens_mdim< vec_t > &out)
Create a conditional probability distribution based on the specified probability distribution.
Random number generator (GSL)
virtual void operator()(vec_t &x) const
Sample the distribution.
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.
double norm
Normalization factor, .
virtual void contour(double level, double theta, vec_t &x)
Return a point on the contour for a specified level given an angle.
rng_gsl rng
Random number generator.
Probability density function based on a histogram.
void set_alt(size_t p_ndim, vec_t &p_peak, mat_t &p_chol, mat_t &p_covar_inv, double p_norm)
Alternate set function for use when covariance matrix has already been decomposed and inverted...
size_t ndim
Number of dimensions.
std::vector< prob_dens_func > list
Vector of one-dimensional distributions.
virtual double pdf(double x) const
The normalized density.
int diagonal_has_zero(const size_t N, mat_t &A)
Return 1 if at least one of the elements in the diagonal is zero.
virtual size_t dim() const
Return the dimensionality.
virtual double operator()() const
Sample from the specified density.
mat_t covar_inv
Inverse of the covariance matrix.
void set_seed(unsigned long int s)
Set the seed.
void set_seed(unsigned long int s)
Set the seed.
std::vector< double > u_low
Lower limits.
static const double x2[5]
size_t ndim
Number of dimensions.
virtual double pdf(double x) const
The normalized density.
virtual double log_pdf(double x) const
The log of the normalized density.
prob_dens_mdim_gaussian()
Create an empty distribution.
const ubvector & partial_sums()
Get reference to partial sums.
double sig_y
The y standard deviation.
Lognormal density function.
virtual double log_pdf(const vec_t &x, const vec_t &x2) const
The log of the conditional probability.
double stddev()
Get the Gaussian width.
mat_t chol
Cholesky decomposition.
void set_center(double cent)
Set the center.
double y0
The y coordinate of the centroid.