[FFmpeg-devel] [PATCH] avutil/tests/lfg.c: added proper normality test
Thomas Turner
thomastdt at googlemail.com
Wed Feb 1 05:56:39 EET 2017
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 <thomastdt at 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
More information about the ffmpeg-devel
mailing list