36 #include <visp3/core/vpConfig.h>
38 #ifdef VISP_HAVE_CATCH2
39 #define CATCH_CONFIG_RUNNER
42 #include <visp3/core/vpGaussRand.h>
43 #include <visp3/core/vpMath.h>
44 #include <visp3/core/vpTime.h>
61 x = a * (x - k * q) - k * r;
72 const long modulo = ntab - 2;
81 for (j = 0; j < ntab; j++) {
93 double ans = (double)y / normalizer;
107 explicit vpUniRandOld(
const long seed = 0)
108 : a(16807), m(2147483647), q(127773), r(2836), normalizer(2147484721.0), x((seed) ? seed : 739806647)
113 virtual ~vpUniRandOld(){};
116 double operator()() {
return draw1(); }
120 TEST_CASE(
"Check Gaussian draw",
"[visp_rand]")
122 std::vector<double> vec(100000);
123 const double sigma = 5.0, mean = -7.5;
128 for (
size_t i = 0; i < vec.size(); i++) {
133 std::cout << vec.size() <<
" Gaussian draw performed in " << chrono.
getDurationMs() <<
" ms" << std::endl;
136 std::cout <<
"Calculated sigma: " << calculated_sigma << std::endl;
137 std::cout <<
"Calculated mean: " << calculated_mean << std::endl;
139 CHECK(calculated_sigma == Approx(sigma).epsilon(0.01));
140 CHECK(calculated_mean == Approx(mean).epsilon(0.01));
143 TEST_CASE(
"Check Gaussian draw independance",
"[visp_rand]")
145 const double sigma = 5.0, mean = -7.5;
147 SECTION(
"Two simultaneous vpGaussRand instances with the same seed should produce the same results")
151 for (
int i = 0; i < 10; i++) {
152 CHECK(rng1() == Approx(rng2()).margin(1e-6));
155 SECTION(
"Two vpGaussRand instances with the same seed should produce the same results")
157 std::vector<double> vec1, vec2;
160 for (
int i = 0; i < 10; i++) {
161 vec1.push_back(rng());
166 for (
int i = 0; i < 10; i++) {
167 vec2.push_back(rng());
170 REQUIRE(vec1.size() == vec2.size());
172 for (
size_t i = 0; i < vec1.size(); i++) {
173 CHECK(vec1[i] == Approx(vec2[i]).margin(1e-6));
178 TEST_CASE(
"Check uniform draw",
"[visp_rand]")
180 const int niters = 500000;
189 for (
int i = 0; i < niters; i++) {
193 if (sqrt(x * x + y * y) <= 1.0) {
197 double pi = 4.0 * inside / niters;
200 double pi_error = pi - M_PI;
201 std::cout <<
"vpUniRand calculated pi: " << pi <<
" in " << chrono.
getDurationMs() <<
" ms" << std::endl;
202 std::cout <<
"pi error: " << pi_error << std::endl;
204 CHECK(pi == Approx(M_PI).margin(0.005));
207 SECTION(
"C++ rand()")
214 for (
int i = 0; i < niters; i++) {
215 double x =
static_cast<double>(rand()) / RAND_MAX;
216 double y =
static_cast<double>(rand()) / RAND_MAX;
218 if (sqrt(x * x + y * y) <= 1.0) {
222 double pi = 4.0 * inside / niters;
225 double pi_error = pi - M_PI;
226 std::cout <<
"C++ rand() calculated pi: " << pi <<
" in " << chrono.
getDurationMs() <<
" ms" << std::endl;
227 std::cout <<
"pi error: " << pi_error << std::endl;
229 CHECK(pi == Approx(M_PI).margin(0.01));
232 SECTION(
"Old ViSP vpUniRand implementation")
239 for (
int i = 0; i < niters; i++) {
243 if (sqrt(x * x + y * y) <= 1.0) {
247 double pi = 4.0 * inside / niters;
250 double pi_error = pi - M_PI;
251 std::cout <<
"Old ViSP vpUniRand implementation calculated pi: " << pi <<
" in " << chrono.
getDurationMs() <<
" ms"
253 std::cout <<
"pi error: " << pi_error << std::endl;
255 CHECK(pi == Approx(M_PI).margin(0.005));
259 TEST_CASE(
"Check uniform draw range",
"[visp_rand]")
261 const int niters = 1000;
264 SECTION(
"Check[0.0, 1.0) range")
266 for (
int i = 0; i < niters; i++) {
273 SECTION(
"Check[-7, 10) range")
275 const int a = -7, b = 10;
276 for (
int i = 0; i < niters; i++) {
283 SECTION(
"Check[-4.5f, 105.3f) range")
285 const float a = -4.5f, b = 105.3f;
286 for (
int i = 0; i < niters; i++) {
293 SECTION(
"Check[14.6, 56.78) range")
295 const double a = 14.6, b = 56.78;
296 for (
int i = 0; i < niters; i++) {
304 TEST_CASE(
"Check uniform draw independance",
"[visp_rand]")
306 SECTION(
"Two simultaneous vpUniRand instances with the same seed should produce the same results")
311 for (
int i = 0; i < 10; i++) {
318 for (
int i = 0; i < 10; i++) {
319 CHECK(rng1.
uniform(-1.0, 1.0) == Approx(rng2.
uniform(-1.0, 1.0)).margin(1e-6));
323 SECTION(
"Two vpUniRand instances with the same seed should produce the same results")
326 std::vector<uint32_t> vec1, vec2;
329 for (
int i = 0; i < 10; i++) {
330 vec1.push_back(rng.
next());
335 for (
int i = 0; i < 10; i++) {
336 vec2.push_back(rng.
next());
339 REQUIRE(vec1.size() == vec2.size());
341 for (
size_t i = 0; i < vec1.size(); i++) {
342 CHECK(vec1[i] == vec2[i]);
346 std::vector<double> vec1, vec2;
349 for (
int i = 0; i < 10; i++) {
350 vec1.push_back(rng.
uniform(-1.0, 2.0));
355 for (
int i = 0; i < 10; i++) {
356 vec2.push_back(rng.
uniform(-1.0, 2.0));
359 REQUIRE(vec1.size() == vec2.size());
361 for (
size_t i = 0; i < vec1.size(); i++) {
362 CHECK(vec1[i] == Approx(vec2[i]).margin(1e-6));
368 int main(
int argc,
char *argv[])
370 Catch::Session session;
373 session.applyCommandLine(argc, argv);
375 int numFailed = session.run();
383 int main() {
return EXIT_SUCCESS; }
void start(bool reset=true)
Class for generating random number with normal probability density.
static double getStdev(const std::vector< double > &v, bool useBesselCorrection=false)
static double getMean(const std::vector< double > &v)
Class for generating random numbers with uniform probability density.
int uniform(int a, int b)