The Chen-Shapiro(CS) test was used to test normality for Lagged Fibonacci PRNG.
Normality Hypothesis Test: The null hypothesis formally tests if the population the sample represents is normally-distributed. For CS, when the normality hypothesis is True, the distribution of QH will have a mean close to 1. Information on CS can be found here: http://www.stata-journal.com/sjpdf.html?articlenum=st0264 http://www.originlab.com/doc/Origin-Help/NormalityTest-Algorithm Signed-off-by: Thomas Turner <thomas...@googlemail.com> --- libavutil/tests/lfg.c | 125 +++++++++++++++++++++++++++++++++++++---------- tests/fate/libavutil.mak | 5 ++ 2 files changed, 104 insertions(+), 26 deletions(-) diff --git a/libavutil/tests/lfg.c b/libavutil/tests/lfg.c index 1425e02..645d39d 100644 --- a/libavutil/tests/lfg.c +++ b/libavutil/tests/lfg.c @@ -20,6 +20,59 @@ #include "libavutil/timer.h" #include "libavutil/lfg.h" +#define TEST 1 + +// Inverse cumulative distribution function +static double inv_cdf(double u) +{ + const double a[4] = { 2.50662823884, + -18.61500062529, + 41.39119773534, + -25.44106049637}; + + const double b[4] = {-8.47351093090, + 23.08336743743, + -21.06224101826, + 3.13082909833}; + + const double c[9] = {0.3374754822726147, + 0.9761690190917186, + 0.1607979714918209, + 0.0276438810333863, + 0.0038405729373609, + 0.0003951896511919, + 0.0000321767881768, + 0.0000002888167364, + 0.0000003960315187}; + + double r; + double x = u - 0.5; + + // Beasley-Springer +#if TEST == 0 + if (fabs(x) < 0.42) { +#endif +#if TEST == 1 + if (u > 0.08 && u < 0.92) { +#endif + double y = x * x; + r = x * (((a[3]*y+a[2])*y+a[1])*y+a[0]) / + ((((b[3]*y+b[2])*y+b[1])*y+b[0])*y+1.0); + } + else {// Moro + r = u; + if (x > 0.0) + r = 1.0 - u; + r = log(-log(r)); + r = c[0] + r*(c[1]+r*(c[2]+r*(c[3]+r*(c[4]+r*(c[5]+r*(c[6]+ + r*(c[7]+r*c[8]))))))); + if (x < 0.0) + r = -r; + } + + return r; +} + int main(void) { int x = 0; @@ -37,38 +90,58 @@ int main(void) } av_log(NULL, AV_LOG_ERROR, "final value:%X\n", x); - /* BMG usage example */ + /* Chen-Shapiro Normality Test */ { - double mean = 1000; - double stddev = 53; - double samp_mean = 0.0, samp_stddev = 0.0; - double samp0, samp1; + double QH = 0.0; + const double stddev = 53.0; + double samp_stddev = 0.0; + const double mean = 1000; + double samp_mean = 0.0; + const int tot_samp = 2000; + double *PRN_arr = av_malloc_array(tot_samp, sizeof(double)); - av_lfg_init(&state, 42); + av_lfg_init(&state, 0x7a69); + if (!PRN_arr) { + fprintf(stderr, "failed to allocate memory!\n"); + return 1; + } - for (i = 0; i < 1000; i += 2) { + for (i = 0; i < tot_samp; i += 2) { double bmg_out[2]; av_bmg_get(&state, bmg_out); - samp0 = bmg_out[0] * stddev + mean; - samp1 = bmg_out[1] * stddev + mean; - samp_mean += samp0 + samp1; - samp_stddev += samp0 * samp0 + samp1 * samp1; - av_log(NULL, AV_LOG_INFO, - "%f\n%f\n", - samp0, - samp1); + PRN_arr[i ] = bmg_out[0] * stddev + mean; + PRN_arr[i+1] = bmg_out[1] * stddev + mean; + samp_mean += PRN_arr[i] + PRN_arr[i+1]; + av_log(NULL, AV_LOG_INFO, "PRN%d : %f\n" + "PRN%d : %f\n", + i, PRN_arr[i], i+1, PRN_arr[i+1]); } - /* TODO: add proper normality test */ - samp_mean /= 1000; - samp_stddev /= 999; - samp_stddev -= (1000.0/999.0)*samp_mean*samp_mean; - samp_stddev = sqrt(samp_stddev); - av_log(NULL, AV_LOG_INFO, "sample mean : %f\n" - "true mean : %f\n" - "sample stddev: %f\n" - "true stddev : %f\n", - samp_mean, mean, samp_stddev, stddev); - } + samp_mean /= tot_samp; + + for (i = 0; i < tot_samp; ++i) { + double temp; + temp = PRN_arr[i] - samp_mean; + temp *= temp; + samp_stddev += temp; + if ( i < (tot_samp - 1)) { + double H_diff; + H_diff = inv_cdf((i + 2.0 - (3.0/8.0)) / (tot_samp + (1.0/4.0))); + H_diff -= inv_cdf((i + 1.0 - (3.0/8.0)) / (tot_samp + (1.0/4.0))); + + QH += ((PRN_arr[i + 1] - PRN_arr[i]) / H_diff); + } + } + samp_stddev = sqrt(samp_stddev / (tot_samp - 1.0)); + QH = 1.0 - QH / ((tot_samp - 1.0) * samp_stddev); + + av_log(NULL, AV_LOG_INFO, "sample mean : %f\n" + "true mean : %f\n" + "sample stddev : %f\n" + "true stddev : %f\n" + "Chen-Shapiro [QH]: %f\n", + samp_mean, mean, samp_stddev, stddev, QH); + av_freep(&PRN_arr); + } return 0; } diff --git a/tests/fate/libavutil.mak b/tests/fate/libavutil.mak index a7bf739..c3b9091 100644 --- a/tests/fate/libavutil.mak +++ b/tests/fate/libavutil.mak @@ -101,6 +101,11 @@ FATE_LIBAVUTIL += fate-imgutils fate-imgutils: libavutil/tests/imgutils$(EXESUF) fate-imgutils: CMD = run libavutil/tests/imgutils +FATE_LIBAVUTIL += fate-lfg +fate-lfg: libavutil/tests/lfg$(EXESUF) +fate-lfg: CMD = run libavutil/tests/lfg +fate-lfg: REF = /dev/null + FATE_LIBAVUTIL += fate-md5 fate-md5: libavutil/tests/md5$(EXESUF) fate-md5: CMD = run libavutil/tests/md5 -- 1.9.1 _______________________________________________ ffmpeg-devel mailing list ffmpeg-devel@ffmpeg.org http://ffmpeg.org/mailman/listinfo/ffmpeg-devel