boost/random/uniform_on_sphere.hpp

   1 /* boost random/uniform_on_sphere.hpp header file
   2  *
   3  * Copyright Jens Maurer 2000-2001
   4  * Copyright Steven Watanabe 2011
   5  * Distributed under the Boost Software License, Version 1.0. (See
   6  * accompanying file LICENSE_1_0.txt or copy at
   7  * http://www.boost.org/LICENSE_1_0.txt)
   8  *
   9  * See http://www.boost.org for most recent version including documentation.
  10  *
  11  * $Id: uniform_on_sphere.hpp 71018 2011-04-05 21:27:52Z steven_watanabe $
  12  *
  13  * Revision history
  14  *  2001-02-18  moved to individual header files
  15  */
  16
  17 #ifndef BOOST_RANDOM_UNIFORM_ON_SPHERE_HPP
  18 #define BOOST_RANDOM_UNIFORM_ON_SPHERE_HPP
  19
  20 #include <vector>
  21 #include <algorithm>     // std::transform
  22 #include <functional>    // std::bind2nd, std::divides
  23 #include <boost/assert.hpp>
  24 #include <boost/random/detail/config.hpp>
  25 #include <boost/random/detail/operators.hpp>
  26 #include <boost/random/normal_distribution.hpp>
  27
  28 namespace boost {
  29 namespace random {
  30
  31 /**
  32  * Instantiations of class template uniform_on_sphere model a
  33  * \random_distribution. Such a distribution produces random
  34  * numbers uniformly distributed on the unit sphere of arbitrary
  35  * dimension @c dim. The @c Cont template parameter must be a STL-like
  36  * container type with begin and end operations returning non-const
  37  * ForwardIterators of type @c Cont::iterator.
  38  */
  39 template<class RealType = double, class Cont = std::vector<RealType> >
  40 class uniform_on_sphere
  41 {
  42 public:
  43     typedef RealType input_type;
  44     typedef Cont result_type;
  45
  46     class param_type
  47     {
  48     public:
  49
  50         typedef uniform_on_sphere distribution_type;
  51
  52         /**
  53          * Constructs the parameters of a uniform_on_sphere
  54          * distribution, given the dimension of the sphere.
  55          */
  56         explicit param_type(int dim_arg = 2) : _dim(dim_arg)
  57         {
  58             BOOST_ASSERT(_dim >= 0);
  59         }
  60
  61         /** Returns the dimension of the sphere. */
  62         int dim() const { return _dim; }
  63
  64         /** Writes the parameters to a @c std::ostream. */
  65         BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, param_type, parm)
  66         {
  67             os << parm._dim;
  68             return os;
  69         }
  70
  71         /** Reads the parameters from a @c std::istream. */
  72         BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, param_type, parm)
  73         {
  74             is >> parm._dim;
  75             return is;
  76         }
  77
  78         /** Returns true if the two sets of parameters are equal. */
  79         BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(param_type, lhs, rhs)
  80         { return lhs._dim == rhs._dim; }
  81
  82         /** Returns true if the two sets of parameters are different. */
  83         BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(param_type)
  84
  85     private:
  86         int _dim;
  87     };
  88
  89     /**
  90      * Constructs a @c uniform_on_sphere distribution.
  91      * @c dim is the dimension of the sphere.
  92      *
  93      * Requires: dim >= 0
  94      */
  95     explicit uniform_on_sphere(int dim_arg = 2)
  96       : _container(dim_arg), _dim(dim_arg) { }
  97
  98     /**
  99      * Constructs a @c uniform_on_sphere distribution from its parameters.
 100      */
 101     explicit uniform_on_sphere(const param_type& parm)
 102       : _container(parm.dim()), _dim(parm.dim()) { }
 103
 104     // compiler-generated copy ctor and assignment operator are fine
 105
 106     /** Returns the dimension of the sphere. */
 107     int dim() const { return _dim; }
 108
 109     /** Returns the parameters of the distribution. */
 110     param_type param() const { return param_type(_dim); }
 111     /** Sets the parameters of the distribution. */
 112     void param(const param_type& parm)
 113     {
 114         _dim = parm.dim();
 115         _container.resize(_dim);
 116     }
 117
 118     /**
 119      * Returns the smallest value that the distribution can produce.
 120      * Note that this is required to approximate the standard library's
 121      * requirements.  The behavior is defined according to lexicographical
 122      * comparison so that for a container type of std::vector,
 123      * dist.min() <= x <= dist.max() where x is any value produced
 124      * by the distribution.
 125      */
 126     result_type min BOOST_PREVENT_MACRO_SUBSTITUTION () const
 127     {
 128         result_type result(_dim);
 129         if(_dim != 0) {
 130             result.front() = RealType(-1.0);
 131         }
 132         return result;
 133     }
 134     /**
 135      * Returns the largest value that the distribution can produce.
 136      * Note that this is required to approximate the standard library's
 137      * requirements.  The behavior is defined according to lexicographical
 138      * comparison so that for a container type of std::vector,
 139      * dist.min() <= x <= dist.max() where x is any value produced
 140      * by the distribution.
 141      */
 142     result_type max BOOST_PREVENT_MACRO_SUBSTITUTION () const
 143     {
 144         result_type result(_dim);
 145         if(_dim != 0) {
 146             result.front() = RealType(1.0);
 147         }
 148         return result;
 149     }
 150
 151     /**
 152      * Effects: Subsequent uses of the distribution do not depend
 153      * on values produced by any engine prior to invoking reset.
 154      */
 155     void reset() { _normal.reset(); }
 156
 157     /**
 158      * Returns a point uniformly distributed over the surface of
 159      * a sphere of dimension dim().
 160      */
 161     template<class Engine>
 162     const result_type & operator()(Engine& eng)
 163     {
 164         RealType sqsum = 0;
 165         for(typename Cont::iterator it = _container.begin();
 166             it != _container.end();
 167             ++it) {
 168             RealType val = _normal(eng);
 169             *it = val;
 170             sqsum += val * val;
 171         }
 172         using std::sqrt;
 173         // for all i: result[i] /= sqrt(sqsum)
 174         std::transform(_container.begin(), _container.end(), _container.begin(),
 175                        std::bind2nd(std::divides<RealType>(), sqrt(sqsum)));
 176         return _container;
 177     }
 178
 179     /**
 180      * Returns a point uniformly distributed over the surface of
 181      * a sphere of dimension param.dim().
 182      */
 183     template<class Engine>
 184     result_type operator()(Engine& eng, const param_type& parm) const
 185     {
 186         return uniform_on_sphere(parm)(eng);
 187     }
 188
 189     /** Writes the distribution to a @c std::ostream. */
 190     BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, uniform_on_sphere, sd)
 191     {
 192         os << sd._dim;
 193         return os;
 194     }
 195
 196     /** Reads the distribution from a @c std::istream. */
 197     BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, uniform_on_sphere, sd)
 198     {
 199         is >> sd._dim;
 200         sd._container.resize(sd._dim);
 201         return is;
 202     }
 203
 204     /**
 205      * Returns true if the two distributions will produce identical
 206      * sequences of values, given equal generators.
 207      */
 208     BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(uniform_on_sphere, lhs, rhs)
 209     { return lhs._dim == rhs._dim && lhs._normal == rhs._normal; }
 210
 211     /**
 212      * Returns true if the two distributions may produce different
 213      * sequences of values, given equal generators.
 214      */
 215     BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(uniform_on_sphere)
 216
 217 private:
 218     normal_distribution<RealType> _normal;
 219     result_type _container;
 220     int _dim;
 221 };
 222
 223 } // namespace random
 224
 225 using random::uniform_on_sphere;
 226
 227 } // namespace boost
 228
 229 #endif // BOOST_RANDOM_UNIFORM_ON_SPHERE_HPP