00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 #include <inttypes.h>
00023 #include <limits.h>
00024 #include <math.h>
00025 #include "lfg.h"
00026 #include "md5.h"
00027 #include "intreadwrite.h"
00028 #include "attributes.h"
00029
00030 void av_cold av_lfg_init(AVLFG *c, unsigned int seed){
00031 uint8_t tmp[16]={0};
00032 int i;
00033
00034 for(i=8; i<64; i+=4){
00035 AV_WL32(tmp, seed); tmp[4]=i;
00036 av_md5_sum(tmp, tmp, 16);
00037 c->state[i ]= AV_RL32(tmp);
00038 c->state[i+1]= AV_RL32(tmp+4);
00039 c->state[i+2]= AV_RL32(tmp+8);
00040 c->state[i+3]= AV_RL32(tmp+12);
00041 }
00042 c->index=0;
00043 }
00044
00045 void av_bmg_get(AVLFG *lfg, double out[2])
00046 {
00047 double x1, x2, w;
00048
00049 do {
00050 x1 = 2.0/UINT_MAX*av_lfg_get(lfg) - 1.0;
00051 x2 = 2.0/UINT_MAX*av_lfg_get(lfg) - 1.0;
00052 w = x1*x1 + x2*x2;
00053 } while (w >= 1.0);
00054
00055 w = sqrt((-2.0 * log(w)) / w);
00056 out[0] = x1 * w;
00057 out[1] = x2 * w;
00058 }
00059
00060 #ifdef TEST
00061 #include "log.h"
00062 #include "timer.h"
00063
00064 int main(void)
00065 {
00066 int x=0;
00067 int i, j;
00068 AVLFG state;
00069
00070 av_lfg_init(&state, 0xdeadbeef);
00071 for (j = 0; j < 10000; j++) {
00072 START_TIMER
00073 for (i = 0; i < 624; i++) {
00074
00075 x+=av_lfg_get(&state);
00076 }
00077 STOP_TIMER("624 calls of av_lfg_get");
00078 }
00079 av_log(NULL, AV_LOG_ERROR, "final value:%X\n", x);
00080
00081
00082 {
00083 double mean = 1000;
00084 double stddev = 53;
00085
00086 av_lfg_init(&state, 42);
00087
00088 for (i = 0; i < 1000; i += 2) {
00089 double bmg_out[2];
00090 av_bmg_get(&state, bmg_out);
00091 av_log(NULL, AV_LOG_INFO,
00092 "%f\n%f\n",
00093 bmg_out[0] * stddev + mean,
00094 bmg_out[1] * stddev + mean);
00095 }
00096 }
00097
00098 return 0;
00099 }
00100 #endif