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 "lfg.h" 24#include "md5.h" 25#include "intreadwrite.h" 26#include "attributes.h" 27 28void av_cold av_lfg_init(AVLFG *c, unsigned int seed){ 29 uint8_t tmp[16]={0}; 30 int i; 31 32 for(i=8; i<64; i+=4){ 33 AV_WL32(tmp, seed); tmp[4]=i; 34 av_md5_sum(tmp, tmp, 16); 35 c->state[i ]= AV_RL32(tmp); 36 c->state[i+1]= AV_RL32(tmp+4); 37 c->state[i+2]= AV_RL32(tmp+8); 38 c->state[i+3]= AV_RL32(tmp+12); 39 } 40 c->index=0; 41} 42 43void av_bmg_get(AVLFG *lfg, double out[2]) 44{ 45 double x1, x2, w; 46 47 do { 48 x1 = 2.0/UINT_MAX*av_lfg_get(lfg) - 1.0; 49 x2 = 2.0/UINT_MAX*av_lfg_get(lfg) - 1.0; 50 w = x1*x1 + x2*x2; 51 } while (w >= 1.0); 52 53 w = sqrt((-2.0 * log(w)) / w); 54 out[0] = x1 * w; 55 out[1] = x2 * w; 56} 57 58#ifdef TEST 59#include "log.h" 60#include "timer.h" 61 62int main(void) 63{ 64 int x=0; 65 int i, j; 66 AVLFG state; 67 68 av_lfg_init(&state, 0xdeadbeef); 69 for (j = 0; j < 10000; j++) { 70 START_TIMER 71 for (i = 0; i < 624; i++) { 72// av_log(NULL,AV_LOG_ERROR, "%X\n", av_lfg_get(&state)); 73 x+=av_lfg_get(&state); 74 } 75 STOP_TIMER("624 calls of av_lfg_get"); 76 } 77 av_log(NULL, AV_LOG_ERROR, "final value:%X\n", x); 78 79 /* BMG usage example */ 80 { 81 double mean = 1000; 82 double stddev = 53; 83 84 av_lfg_init(&state, 42); 85 86 for (i = 0; i < 1000; i += 2) { 87 double bmg_out[2]; 88 av_bmg_get(&state, bmg_out); 89 av_log(NULL, AV_LOG_INFO, 90 "%f\n%f\n", 91 bmg_out[0] * stddev + mean, 92 bmg_out[1] * stddev + mean); 93 } 94 } 95 96 return 0; 97} 98#endif 99