// Copyright (C) 2007 Davis E. King (davis@dlib.net)
// License: Boost Software License See LICENSE.txt for the full license.
#ifndef DLIB_RAND_KERNEl_1_
#define DLIB_RAND_KERNEl_1_
#include <string>
#include <complex>
#include "../algs.h"
#include "rand_kernel_abstract.h"
#include "mersenne_twister.h"
#include "../is_kind.h"
#include <iostream>
#include "../serialize.h"
#include "../string.h"
namespace dlib
{
class rand
{
/*!
INITIAL VALUE
- seed == ""
CONVENTION
- the random numbers come from the boost mersenne_twister code
- get_seed() == seed
!*/
public:
// These typedefs are here for backwards compatibility with older versions of dlib.
typedef rand kernel_1a;
typedef rand float_1a;
rand(
)
{
init();
}
rand (
time_t seed_value
)
{
init();
set_seed(cast_to_string(seed_value));
}
rand (
const std::string& seed_value
)
{
init();
set_seed(seed_value);
}
virtual ~rand(
)
{}
void clear(
)
{
mt.seed();
seed.clear();
has_gaussian = false;
next_gaussian = 0;
// prime the generator a bit
for (int i = 0; i < 10000; ++i)
mt();
}
const std::string& get_seed (
)
{
return seed;
}
void set_seed (
const std::string& value
)
{
seed = value;
// make sure we do the seeding so that using a seed of "" gives the same
// state as calling this->clear()
if (value.size() != 0)
{
uint32 s = 0;
for (std::string::size_type i = 0; i < seed.size(); ++i)
{
s = (s*37) + static_cast<uint32>(seed[i]);
}
mt.seed(s);
}
else
{
mt.seed();
}
// prime the generator a bit
for (int i = 0; i < 10000; ++i)
mt();
has_gaussian = false;
next_gaussian = 0;
}
unsigned char get_random_8bit_number (
)
{
return static_cast<unsigned char>(mt());
}
uint16 get_random_16bit_number (
)
{
return static_cast<uint16>(mt());
}
inline uint32 get_random_32bit_number (
)
{
return mt();
}
inline uint64 get_random_64bit_number (
)
{
const uint64 a = get_random_32bit_number();
const uint64 b = get_random_32bit_number();
return (a<<32)|b;
}
double get_double_in_range (
double begin,
double end
)
{
DLIB_ASSERT(begin <= end);
return begin + get_random_double()*(end-begin);
}
long long get_integer_in_range(
long long begin,
long long end
)
{
DLIB_ASSERT(begin <= end);
if (begin == end)
return begin;
auto r = get_random_64bit_number();
const auto limit = std::numeric_limits<decltype(r)>::max();
const auto range = end-begin;
// Use rejection sampling to remove the biased sampling you would get with
// the naive get_random_64bit_number()%range sampling.
while(r >= (limit/range)*range)
r = get_random_64bit_number();
return begin + static_cast<long long>(r%range);
}
long long get_integer(
long long end
)
{
DLIB_ASSERT(end >= 0);
return get_integer_in_range(0,end);
}
double get_random_double (
)
{
uint32 temp;
temp = rand::get_random_32bit_number();
temp &= 0xFFFFFF;
double val = static_cast<double>(temp);
val *= 0x1000000;
temp = rand::get_random_32bit_number();
temp &= 0xFFFFFF;
val += temp;
val /= max_val;
if (val < 1.0)
{
return val;
}
else
{
// return a value slightly less than 1.0
return 1.0 - std::numeric_limits<double>::epsilon();
}
}
float get_random_float (
)
{
uint32 temp;
temp = rand::get_random_32bit_number();
temp &= 0xFFFFFF;
const float scale = 1.0/0x1000000;
const float val = static_cast<float>(temp)*scale;
if (val < 1.0f)
{
return val;
}
else
{
// return a value slightly less than 1.0
return 1.0f - std::numeric_limits<float>::epsilon();
}
}
std::complex<double> get_random_complex_gaussian (
)
{
double x1, x2, w;
const double rndmax = std::numeric_limits<dlib::uint32>::max();
// Generate a pair of Gaussian random numbers using the Box-Muller transformation.
do
{
const double rnd1 = get_random_32bit_number()/rndmax;
const double rnd2 = get_random_32bit_number()/rndmax;
x1 = 2.0 * rnd1 - 1.0;
x2 = 2.0 * rnd2 - 1.0;
w = x1 * x1 + x2 * x2;
} while ( w >= 1.0 );
w = std::sqrt( (-2.0 * std::log( w ) ) / w );
return std::complex<double>(x1 * w, x2 * w);
}
double get_random_gaussian (
)
{
if (has_gaussian)
{
has_gaussian = false;
return next_gaussian;
}
std::complex<double> r = get_random_complex_gaussian();
next_gaussian = r.imag();
has_gaussian = true;
return r.real();
}
double get_random_exponential (
double lambda
)
{
DLIB_ASSERT(lambda > 0, "lambda must be greater than zero");
double u = 0.0;
while (u == 0.0)
u = get_random_double();
return -std::log( u ) / lambda;
}
double get_random_weibull (
double lambda,
double k,
double gamma
)
{
DLIB_ASSERT(k > 0, "k must be greater than zero");
DLIB_ASSERT(lambda > 0, "lambda must be greater than zero");
double u = 0.0;
while (u == 0.0)
u = get_random_double();
return gamma + lambda*std::pow(-std::log(u), 1.0 / k);
}
double get_random_beta (
double alpha,
double beta
)
{
DLIB_CASSERT(alpha > 0, "alpha must be greater than zero")
DLIB_CASSERT(beta > 0, "beta must be greater than zero");
auto u = std::pow(get_random_double(), 1 / alpha);
auto v = std::pow(get_random_double(), 1 / beta);
while ((u + v) > 1 || (u == 0 && v == 0))
{
u = std::pow(get_random_double(), 1 / alpha);
v = std::pow(get_random_double(), 1 / beta);
}
return u / (u + v);
}
void swap (
rand& item
)
{
exchange(mt,item.mt);
exchange(seed, item.seed);
exchange(has_gaussian, item.has_gaussian);
exchange(next_gaussian, item.next_gaussian);
}
friend void serialize(
const rand& item,
std::ostream& out
);
friend void deserialize(
rand& item,
std::istream& in
);
private:
void init()
{
// prime the generator a bit
for (int i = 0; i < 10000; ++i)
mt();
max_val = 0xFFFFFF;
max_val *= 0x1000000;
max_val += 0xFFFFFF;
max_val += 0.05;
has_gaussian = false;
next_gaussian = 0;
}
mt19937 mt;
std::string seed;
double max_val;
bool has_gaussian;
double next_gaussian;
};
inline void swap (
rand& a,
rand& b
) { a.swap(b); }
template <>
struct is_rand<rand>
{
static const bool value = true;
};
inline void serialize(
const rand& item,
std::ostream& out
)
{
int version = 1;
serialize(version, out);
serialize(item.mt, out);
serialize(item.seed, out);
serialize(item.has_gaussian, out);
serialize(item.next_gaussian, out);
}
inline void deserialize(
rand& item,
std::istream& in
)
{
int version;
deserialize(version, in);
if (version != 1)
throw serialization_error("Error deserializing object of type rand: unexpected version.");
deserialize(item.mt, in);
deserialize(item.seed, in);
deserialize(item.has_gaussian, in);
deserialize(item.next_gaussian, in);
}
}
#endif // DLIB_RAND_KERNEl_1_