144 lines
4.4 KiB
C++
144 lines
4.4 KiB
C++
// SPDX-License-Identifier: BSD-2-Clause
|
|
|
|
// This code is part of the sfizz library and is licensed under a BSD 2-clause
|
|
// license. You should have receive a LICENSE.md file along with the code.
|
|
// If not, contact the sfizz maintainers at https://github.com/sfztools/sfizz
|
|
|
|
#include "catch2/catch.hpp"
|
|
#include "sfizz/MathHelpers.h"
|
|
#include <vector>
|
|
|
|
template <class T>
|
|
static inline T squared(T x)
|
|
{
|
|
return x * x;
|
|
}
|
|
|
|
/**
|
|
* @brief Check the behavior of a uniform real random generator
|
|
*
|
|
* - ensure all results to be in range [min;max]
|
|
* - ensure to have at least one element in every division of the result range
|
|
*/
|
|
template <class Rand, class Dist>
|
|
static bool uniformRandomTest(float min, float max, unsigned gen, unsigned div)
|
|
{
|
|
fast_rand prng;
|
|
fast_real_distribution<float> dist(min, max);
|
|
std::vector<unsigned> counts(div);
|
|
|
|
for (unsigned i = 0; i < gen; ++i) {
|
|
float r = dist(prng);
|
|
if (r < min || r > max)
|
|
return false;
|
|
unsigned d = clamp<int>(div * (r - min) / (max - min), 0, div - 1);
|
|
++counts[d];
|
|
}
|
|
|
|
for (unsigned count : counts) {
|
|
if (count == 0)
|
|
return false;
|
|
}
|
|
return true;
|
|
}
|
|
|
|
TEST_CASE("[Random] Fast random generation")
|
|
{
|
|
unsigned numGenerations = 1000;
|
|
unsigned numDivisions = 128;
|
|
|
|
auto& test = uniformRandomTest<fast_rand, fast_real_distribution<float>>;
|
|
|
|
REQUIRE(test(0.0f, 1.0f, numGenerations, numDivisions));
|
|
REQUIRE(test(-1.0f, 1.0f, numGenerations, numDivisions));
|
|
REQUIRE(test(0.0f, 123.0f, numGenerations, numDivisions));
|
|
REQUIRE(test(-123.0f, 0.0f, numGenerations, numDivisions));
|
|
}
|
|
|
|
/**
|
|
* @brief Check the behavior of a gaussian real random generator
|
|
*
|
|
* - compute a histogram of random generations
|
|
* - compare it against the standard library normal distribution
|
|
*/
|
|
template <unsigned Quality>
|
|
static bool gaussianRandomTest(double mean, float variance, size_t numGen, size_t histSize, double maxAbsErr)
|
|
{
|
|
fast_gaussian_generator<float, Quality> gen(mean, variance);
|
|
|
|
std::mt19937_64 prng;
|
|
std::normal_distribution<float> dist(mean, variance);
|
|
|
|
// the tested bounds
|
|
const double min = -1.0 + mean;
|
|
const double max = +1.0 + mean;
|
|
|
|
// generate, quantify, count occurrences
|
|
std::vector<size_t> counts(histSize);
|
|
for (size_t i = 0; i < numGen; ++i) {
|
|
double value = gen();
|
|
double normalized = (value - min) / (max - min);
|
|
long bin = std::lround(histSize * normalized);
|
|
if (bin >= 0 && static_cast<size_t>(bin) < histSize)
|
|
counts[bin] += 1;
|
|
}
|
|
std::vector<size_t> referenceCounts(histSize);
|
|
for (size_t i = 0; i < numGen; ++i) {
|
|
double value = dist(prng);
|
|
double normalized = (value - min) / (max - min);
|
|
long bin = std::lround(histSize * normalized);
|
|
if (bin >= 0 && static_cast<size_t>(bin) < histSize)
|
|
referenceCounts[bin] += 1;
|
|
}
|
|
|
|
// bin probabilities
|
|
std::vector<double> proba(histSize);
|
|
std::vector<double> referenceProba(histSize);
|
|
for (size_t i = 0; i < histSize; ++i) {
|
|
proba[i] = static_cast<double>(counts[i]) / numGen;
|
|
referenceProba[i] = static_cast<double>(referenceCounts[i]) / numGen;
|
|
}
|
|
|
|
// normalize to unity, for the sake of comparison
|
|
{
|
|
double sum = 0.0;
|
|
for (double x : referenceProba)
|
|
sum += x;
|
|
double scaleFactor = 1.0 / sum;
|
|
for (size_t i = 0; i < histSize; ++i) {
|
|
proba[i] *= scaleFactor;
|
|
referenceProba[i] *= scaleFactor;
|
|
}
|
|
}
|
|
|
|
auto normalizeToUnity = [](absl::Span<double> v) {
|
|
double s = 0;
|
|
for (double x : v)
|
|
s += x;
|
|
for (double& x : v)
|
|
x *= 1.0 / s;
|
|
};
|
|
normalizeToUnity(absl::MakeSpan(proba));
|
|
normalizeToUnity(absl::MakeSpan(referenceProba));
|
|
|
|
// compare
|
|
for (size_t i = 0; i < histSize; ++i) {
|
|
double absErr = std::abs(proba[i] - referenceProba[i]);
|
|
if (absErr > maxAbsErr)
|
|
return false;
|
|
}
|
|
|
|
return true;
|
|
}
|
|
|
|
TEST_CASE("[Random] Gaussian random generation")
|
|
{
|
|
unsigned numGenerations = 16384;
|
|
unsigned numDivisions = 128;
|
|
|
|
const double maxAbsErr = 0.05; // PDF ±5%
|
|
|
|
REQUIRE(gaussianRandomTest<4>(0.0f, 0.25f, numGenerations, numDivisions, maxAbsErr));
|
|
REQUIRE(gaussianRandomTest<4>(0.0f, 0.50f, numGenerations, numDivisions, maxAbsErr));
|
|
REQUIRE(gaussianRandomTest<4>(0.0f, 0.75f, numGenerations, numDivisions, maxAbsErr));
|
|
}
|