Subversion Repositories Kolibri OS

Rev

Blame | Last modification | View Log | RSS feed

  1. /*
  2.  * Lagged Fibonacci PRNG
  3.  * Copyright (c) 2008 Michael Niedermayer
  4.  *
  5.  * This file is part of FFmpeg.
  6.  *
  7.  * FFmpeg is free software; you can redistribute it and/or
  8.  * modify it under the terms of the GNU Lesser General Public
  9.  * License as published by the Free Software Foundation; either
  10.  * version 2.1 of the License, or (at your option) any later version.
  11.  *
  12.  * FFmpeg is distributed in the hope that it will be useful,
  13.  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  14.  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  15.  * Lesser General Public License for more details.
  16.  *
  17.  * You should have received a copy of the GNU Lesser General Public
  18.  * License along with FFmpeg; if not, write to the Free Software
  19.  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
  20.  */
  21.  
  22. #include <inttypes.h>
  23. #include <limits.h>
  24. #include <math.h>
  25. #include "lfg.h"
  26. #include "md5.h"
  27. #include "intreadwrite.h"
  28. #include "attributes.h"
  29.  
  30. av_cold void av_lfg_init(AVLFG *c, unsigned int seed)
  31. {
  32.     uint8_t tmp[16] = { 0 };
  33.     int i;
  34.  
  35.     for (i = 8; i < 64; i += 4) {
  36.         AV_WL32(tmp, seed);
  37.         tmp[4] = i;
  38.         av_md5_sum(tmp, tmp, 16);
  39.         c->state[i    ] = AV_RL32(tmp);
  40.         c->state[i + 1] = AV_RL32(tmp + 4);
  41.         c->state[i + 2] = AV_RL32(tmp + 8);
  42.         c->state[i + 3] = AV_RL32(tmp + 12);
  43.     }
  44.     c->index = 0;
  45. }
  46.  
  47. void av_bmg_get(AVLFG *lfg, double out[2])
  48. {
  49.     double x1, x2, w;
  50.  
  51.     do {
  52.         x1 = 2.0 / UINT_MAX * av_lfg_get(lfg) - 1.0;
  53.         x2 = 2.0 / UINT_MAX * av_lfg_get(lfg) - 1.0;
  54.         w  = x1 * x1 + x2 * x2;
  55.     } while (w >= 1.0);
  56.  
  57.     w = sqrt((-2.0 * log(w)) / w);
  58.     out[0] = x1 * w;
  59.     out[1] = x2 * w;
  60. }
  61.  
  62. #ifdef TEST
  63. #include "log.h"
  64. #include "timer.h"
  65.  
  66. int main(void)
  67. {
  68.     int x = 0;
  69.     int i, j;
  70.     AVLFG state;
  71.  
  72.     av_lfg_init(&state, 0xdeadbeef);
  73.     for (j = 0; j < 10000; j++) {
  74.         START_TIMER
  75.         for (i = 0; i < 624; i++) {
  76.             //av_log(NULL, AV_LOG_ERROR, "%X\n", av_lfg_get(&state));
  77.             x += av_lfg_get(&state);
  78.         }
  79.         STOP_TIMER("624 calls of av_lfg_get");
  80.     }
  81.     av_log(NULL, AV_LOG_ERROR, "final value:%X\n", x);
  82.  
  83.     /* BMG usage example */
  84.     {
  85.         double mean   = 1000;
  86.         double stddev = 53;
  87.  
  88.         av_lfg_init(&state, 42);
  89.  
  90.         for (i = 0; i < 1000; i += 2) {
  91.             double bmg_out[2];
  92.             av_bmg_get(&state, bmg_out);
  93.             av_log(NULL, AV_LOG_INFO,
  94.                    "%f\n%f\n",
  95.                    bmg_out[0] * stddev + mean,
  96.                    bmg_out[1] * stddev + mean);
  97.         }
  98.     }
  99.  
  100.     return 0;
  101. }
  102. #endif
  103.