00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #ifndef BZ_RANDOM_NORMAL
00022 #define BZ_RANDOM_NORMAL
00023
00024 #ifndef BZ_RANDOM_UNIFORM
00025 #include <random/uniform.h>
00026 #endif
00027
00028 BZ_NAMESPACE(ranlib)
00029
00030 template<typename T = double, typename IRNG = defaultIRNG,
00031 typename stateTag = defaultState>
00032 class NormalUnit : public UniformOpen<T,IRNG,stateTag>
00033 {
00034 public:
00035 typedef T T_numtype;
00036
00037 T random()
00038 {
00039 const T s = 0.449871, t = -0.386595, a = 0.19600, b = 0.25472;
00040 const T r1 = 0.27597, r2 = 0.27846;
00041
00042 T u, v;
00043
00044 for (;;) {
00045
00046
00047
00048
00049
00050
00051 u = this->getUniform();
00052 v = 1.715527769921413592960379282557544956242L
00053 * (this->getUniform() - 0.5);
00054
00055
00056 T x = u - s;
00057 T y = fabs(v) - t;
00058 T q = x*x + y*(a*y - b*x);
00059
00060
00061 if (q < r1)
00062 break;
00063
00064
00065 if (q > r2)
00066 continue;
00067
00068
00069 if (v*v <= -4.0 * log(u)*u*u)
00070 break;
00071 }
00072
00073 return v/u;
00074 }
00075
00076 };
00077
00078
00079 template<typename T = double, typename IRNG = defaultIRNG,
00080 typename stateTag = defaultState>
00081 class Normal : public NormalUnit<T,IRNG,stateTag> {
00082
00083 public:
00084 typedef T T_numtype;
00085
00086 Normal(T mean, T standardDeviation)
00087 {
00088 mean_ = mean;
00089 standardDeviation_ = standardDeviation;
00090 }
00091
00092 T random()
00093 {
00094 return mean_ + standardDeviation_
00095 * NormalUnit<T,IRNG,stateTag>::random();
00096 }
00097
00098 private:
00099 T mean_;
00100 T standardDeviation_;
00101 };
00102
00103 BZ_NAMESPACE_END
00104
00105 #endif // BZ_RANDOM_NORMAL