/usr/include/boost/random
Edit: /usr/include/boost/random/uniform_on_sphere.hpp (8812B)
/* boost random/uniform_on_sphere.hpp header file
*
* Copyright Jens Maurer 2000-2001
* Copyright Steven Watanabe 2011
* Distributed under the Boost Software License, Version 1.0. (See
* accompanying file LICENSE_1_0.txt or copy at
* http://www.boost.org/LICENSE_1_0.txt)
*
* See http://www.boost.org for most recent version including documentation.
*
* $Id$
*
* Revision history
* 2001-02-18 moved to individual header files
*/
#ifndef BOOST_RANDOM_UNIFORM_ON_SPHERE_HPP
#define BOOST_RANDOM_UNIFORM_ON_SPHERE_HPP
#include
#include // std::transform
#include // std::bind2nd, std::divides
#include
#include
#include
#include
namespace boost {
namespace random {
/**
* Instantiations of class template uniform_on_sphere model a
* \random_distribution. Such a distribution produces random
* numbers uniformly distributed on the unit sphere of arbitrary
* dimension @c dim. The @c Cont template parameter must be a STL-like
* container type with begin and end operations returning non-const
* ForwardIterators of type @c Cont::iterator.
*/
template >
class uniform_on_sphere
{
public:
typedef RealType input_type;
typedef Cont result_type;
class param_type
{
public:
typedef uniform_on_sphere distribution_type;
/**
* Constructs the parameters of a uniform_on_sphere
* distribution, given the dimension of the sphere.
*/
explicit param_type(int dim_arg = 2) : _dim(dim_arg)
{
BOOST_ASSERT(_dim >= 0);
}
/** Returns the dimension of the sphere. */
int dim() const { return _dim; }
/** Writes the parameters to a @c std::ostream. */
BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, param_type, parm)
{
os << parm._dim;
return os;
}
/** Reads the parameters from a @c std::istream. */
BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, param_type, parm)
{
is >> parm._dim;
return is;
}
/** Returns true if the two sets of parameters are equal. */
BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(param_type, lhs, rhs)
{ return lhs._dim == rhs._dim; }
/** Returns true if the two sets of parameters are different. */
BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(param_type)
private:
int _dim;
};
/**
* Constructs a @c uniform_on_sphere distribution.
* @c dim is the dimension of the sphere.
*
* Requires: dim >= 0
*/
explicit uniform_on_sphere(int dim_arg = 2)
: _container(dim_arg), _dim(dim_arg) { }
/**
* Constructs a @c uniform_on_sphere distribution from its parameters.
*/
explicit uniform_on_sphere(const param_type& parm)
: _container(parm.dim()), _dim(parm.dim()) { }
// compiler-generated copy ctor and assignment operator are fine
/** Returns the dimension of the sphere. */
int dim() const { return _dim; }
/** Returns the parameters of the distribution. */
param_type param() const { return param_type(_dim); }
/** Sets the parameters of the distribution. */
void param(const param_type& parm)
{
_dim = parm.dim();
_container.resize(_dim);
}
/**
* Returns the smallest value that the distribution can produce.
* Note that this is required to approximate the standard library's
* requirements. The behavior is defined according to lexicographical
* comparison so that for a container type of std::vector,
* dist.min() <= x <= dist.max() where x is any value produced
* by the distribution.
*/
result_type min BOOST_PREVENT_MACRO_SUBSTITUTION () const
{
result_type result(_dim);
if(_dim != 0) {
result.front() = RealType(-1.0);
}
return result;
}
/**
* Returns the largest value that the distribution can produce.
* Note that this is required to approximate the standard library's
* requirements. The behavior is defined according to lexicographical
* comparison so that for a container type of std::vector,
* dist.min() <= x <= dist.max() where x is any value produced
* by the distribution.
*/
result_type max BOOST_PREVENT_MACRO_SUBSTITUTION () const
{
result_type result(_dim);
if(_dim != 0) {
result.front() = RealType(1.0);
}
return result;
}
/**
* Effects: Subsequent uses of the distribution do not depend
* on values produced by any engine prior to invoking reset.
*/
void reset() {}
/**
* Returns a point uniformly distributed over the surface of
* a sphere of dimension dim().
*/
template
const result_type & operator()(Engine& eng)
{
using std::sqrt;
switch(_dim)
{
case 0: break;
case 1:
{
if(uniform_01()(eng) < 0.5) {
*_container.begin() = -1;
} else {
*_container.begin() = 1;
}
break;
}
case 2:
{
uniform_01 uniform;
RealType sqsum;
RealType x, y;
do {
x = uniform(eng) * 2 - 1;
y = uniform(eng) * 2 - 1;
sqsum = x*x + y*y;
} while(sqsum == 0 || sqsum > 1);
RealType mult = 1/sqrt(sqsum);
typename Cont::iterator iter = _container.begin();
*iter = x * mult;
iter++;
*iter = y * mult;
break;
}
case 3:
{
uniform_01 uniform;
RealType sqsum;
RealType x, y;
do {
x = uniform(eng) * 2 - 1;
y = uniform(eng) * 2 - 1;
sqsum = x*x + y*y;
} while(sqsum > 1);
RealType mult = 2 * sqrt(1 - sqsum);
typename Cont::iterator iter = _container.begin();
*iter = x * mult;
++iter;
*iter = y * mult;
++iter;
*iter = 2 * sqsum - 1;
break;
}
default:
{
detail::unit_normal_distribution normal;
RealType sqsum;
do {
sqsum = 0;
for(typename Cont::iterator it = _container.begin();
it != _container.end();
++it) {
RealType val = normal(eng);
*it = val;
sqsum += val * val;
}
} while(sqsum == 0);
// for all i: result[i] /= sqrt(sqsum)
RealType inverse_distance = 1 / sqrt(sqsum);
for(typename Cont::iterator it = _container.begin();
it != _container.end();
++it) {
*it *= inverse_distance;
}
}
}
return _container;
}
/**
* Returns a point uniformly distributed over the surface of
* a sphere of dimension param.dim().
*/
template
result_type operator()(Engine& eng, const param_type& parm) const
{
return uniform_on_sphere(parm)(eng);
}
/** Writes the distribution to a @c std::ostream. */
BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, uniform_on_sphere, sd)
{
os << sd._dim;
return os;
}
/** Reads the distribution from a @c std::istream. */
BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, uniform_on_sphere, sd)
{
is >> sd._dim;
sd._container.resize(sd._dim);
return is;
}
/**
* Returns true if the two distributions will produce identical
* sequences of values, given equal generators.
*/
BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(uniform_on_sphere, lhs, rhs)
{ return lhs._dim == rhs._dim; }
/**
* Returns true if the two distributions may produce different
* sequences of values, given equal generators.
*/
BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(uniform_on_sphere)
private:
result_type _container;
int _dim;
};
} // namespace random
using random::uniform_on_sphere;
} // namespace boost
#endif // BOOST_RANDOM_UNIFORM_ON_SPHERE_HPP