123 lines
2.8 KiB
C++
123 lines
2.8 KiB
C++
#ifndef RANDOM_NUMBER_GENERATORS_HPP
|
|
#define RANDOM_NUMBER_GENERATORS_HPP
|
|
#include <cmath>
|
|
#include <vector>
|
|
#include <cassert>
|
|
#include <iostream>
|
|
|
|
#define WEIGHTED_RANDOM_DEBUG 0
|
|
|
|
namespace RandomNumberGenerators {
|
|
|
|
inline float uniform()
|
|
/* Uniform random number generator x(n+1)= a*x(n) mod c
|
|
with a = pow(7,5) and c = pow(2,31)-1.
|
|
Copyright (c) Tao Pang 1997. */
|
|
{
|
|
const int ia=16807,ic=2147483647,iq=127773,ir=2836;
|
|
int il,ih,it;
|
|
float rc;
|
|
static int iseed = rand();
|
|
ih = iseed/iq;
|
|
il = iseed%iq;
|
|
it = ia*il-ir*ih;
|
|
if (it > 0)
|
|
{
|
|
iseed = it;
|
|
}
|
|
else
|
|
{
|
|
iseed = ic+it;
|
|
}
|
|
rc = ic;
|
|
return iseed/rc;
|
|
}
|
|
|
|
inline float gaussian(float mean, float sigma)
|
|
{
|
|
|
|
float x1, x2, w, y1, y2;
|
|
|
|
do {
|
|
x1 = 2.0 * uniform() - 1.0;
|
|
x2 = 2.0 * uniform() - 1.0;
|
|
w = x1 * x1 + x2 * x2;
|
|
} while ( w >= 1.0 );
|
|
|
|
w = sqrt( (-2.0 * log( w ) ) / w );
|
|
y1 = x1 * w;
|
|
y2 = x2 * w;
|
|
|
|
float ret = y1*sigma + mean;
|
|
|
|
return ret;
|
|
}
|
|
|
|
inline std::size_t uniformInteger(std::size_t upperBound=1) {
|
|
|
|
|
|
/// @bug there was a man entry about how this leads to a lousy uniform
|
|
/// @bug distribution in practice. should probably review
|
|
assert(upperBound > 0);
|
|
return ((rand()) % ((int)upperBound));
|
|
}
|
|
|
|
|
|
|
|
|
|
/// Randomizes from probabilistically weighted distribution. Thus,
|
|
/// sum of passed in weights should be 1.0
|
|
inline std::size_t weightedRandomNormalized(std::vector<float> weights) {
|
|
|
|
// Choose a random bounded mass between 0 and 1
|
|
float cutoff = ((float)(rand())) / (float)RAND_MAX;
|
|
|
|
//std::cout << "cutoff : " << cutoff << std::endl;
|
|
|
|
// Sum up mass, stopping when cutoff is reached. This is the typical
|
|
// weighted sampling algorithm.
|
|
float mass = 0;
|
|
for (std::size_t i = 0; i< weights.size() ; i++) {
|
|
mass += weights[i];
|
|
//std::cout << "mass: " << mass << std::endl;
|
|
if (mass >= cutoff)
|
|
return i;
|
|
}
|
|
|
|
// Just in case something slips through the cracks
|
|
return weights.size()-1;
|
|
}
|
|
|
|
inline std::size_t weightedRandom(const std::vector<int> & weights, unsigned int weightTotalHint = 0) {
|
|
|
|
|
|
if (weightTotalHint == 0) {
|
|
for (std::size_t i = 0; i < weights.size();i++)
|
|
weightTotalHint += weights[i];
|
|
}
|
|
|
|
const int sampledSum = uniformInteger(weightTotalHint);
|
|
int sum = 0;
|
|
if (WEIGHTED_RANDOM_DEBUG) std::cout << "[RNG::weightedRandom()] weightTotal = " << weightTotalHint <<
|
|
std::endl;
|
|
|
|
for (std::size_t i = 0; i < weights.size();i++) {
|
|
if (WEIGHTED_RANDOM_DEBUG)
|
|
std::cout << "[RNG::weightedRandom()] weight[" << i << "] = " << weights[i] <<
|
|
std::endl;
|
|
|
|
sum += weights[i];
|
|
if (sampledSum <= sum) {
|
|
if (WEIGHTED_RANDOM_DEBUG)
|
|
std::cout << "[RNG::weightedRandom()] sampled index " << i << "(" <<
|
|
"running sum = " << sum << ", sampled sum = " << sampledSum << std::endl;
|
|
return i;
|
|
}
|
|
}
|
|
|
|
return weights.size()-1;
|
|
}
|
|
|
|
}
|
|
#endif
|