From patchwork Wed Feb 1 03:56:39 2017 Content-Type: text/plain; charset="utf-8" MIME-Version: 1.0 Content-Transfer-Encoding: 7bit X-Patchwork-Submitter: Thomas Turner X-Patchwork-Id: 2394 Delivered-To: ffmpegpatchwork@gmail.com Received: by 10.103.89.21 with SMTP id n21csp2258396vsb; Tue, 31 Jan 2017 20:04:27 -0800 (PST) X-Received: by 10.28.129.147 with SMTP id c141mr23018577wmd.12.1485921867802; Tue, 31 Jan 2017 20:04:27 -0800 (PST) Return-Path: Received: from ffbox0-bg.mplayerhq.hu (ffbox0-bg.ffmpeg.org. [79.124.17.100]) by mx.google.com with ESMTP id b72si19711735wmf.22.2017.01.31.20.04.25; Tue, 31 Jan 2017 20:04:27 -0800 (PST) Received-SPF: pass (google.com: domain of ffmpeg-devel-bounces@ffmpeg.org designates 79.124.17.100 as permitted sender) client-ip=79.124.17.100; Authentication-Results: mx.google.com; dkim=neutral (body hash did not verify) header.i=@googlemail.com; spf=pass (google.com: domain of ffmpeg-devel-bounces@ffmpeg.org designates 79.124.17.100 as permitted sender) smtp.mailfrom=ffmpeg-devel-bounces@ffmpeg.org; dmarc=fail (p=QUARANTINE sp=QUARANTINE dis=NONE) header.from=googlemail.com Received: from [127.0.1.1] (localhost [127.0.0.1]) by ffbox0-bg.mplayerhq.hu (Postfix) with ESMTP id 8DEF86891F1; Wed, 1 Feb 2017 06:04:20 +0200 (EET) X-Original-To: ffmpeg-devel@ffmpeg.org Delivered-To: ffmpeg-devel@ffmpeg.org Received: from mail-oi0-f65.google.com (mail-oi0-f65.google.com [209.85.218.65]) by ffbox0-bg.mplayerhq.hu (Postfix) with ESMTPS id 8C193680D23 for ; Wed, 1 Feb 2017 06:04:14 +0200 (EET) Received: by mail-oi0-f65.google.com with SMTP id w144so30021949oiw.1 for ; Tue, 31 Jan 2017 20:04:17 -0800 (PST) DKIM-Signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=googlemail.com; s=20161025; h=from:to:subject:date:message-id; bh=HqljGL/pwXWPzTiP+Swn8J1WJ25eIFz+EMEj03ry20k=; b=NWdWLGrc5xGHrArpbcwUX+YOL6M9pVjTcszf5htRT8wL0Nxp9TG/FcD/davvIE3Oya Ac55ZqyYj+Epsjb0nzreniXwa2qCXXOq29u4zpNZ/VzwLOh54YJgXNlRamgjntuRW5mE Mf5Tds8gH2z0DRpbL3C8nOaMlvmRuEH7cUfiohZTljuYYEiBjacUQPj6xVwGZq+APHYA A2rCZOPuAbywdney5aKPeUtDBP/lBOpVCWy/8E5XZlh8Eo7YgA9+91Vc34iSP3xZ/oV6 ACBn5UCzwsiQohG4cm4TLLOk4Lmp8SwFyzlUJSkia93Qcb0M5jpS+JGnRvxU7NwLb1wl wtRA== X-Google-DKIM-Signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=1e100.net; s=20161025; h=x-gm-message-state:from:to:subject:date:message-id; bh=HqljGL/pwXWPzTiP+Swn8J1WJ25eIFz+EMEj03ry20k=; b=YDhjmNfRU/ofD0lfgj3uMslYXrPvA/JWVJn5reKeVVq+yBcjSv3bMp4d4P6KarLm28 zLB3TW7hb0x8Bc1aIGt/kxNECKAScdT/STjESOOXWfhfjpTDxBIlxvtPYDNRK2kUwKfz 8X+gfQUqbL8pyYmd3WdbuAUoYkjqCeZHjAFk1Wj/h8/w5pADkBVqNxGD6aulxfg6C7KN HFYgFSpWp+9VtlBmeu8Z/C88P0Z1WRl0a0RDjRZ7uUehhysT/h5Bm5QJniv76XCzkTVb idLqNFtgUQxjivl5UlwAF5iC3IMxmmvzdjjGV0a4TJHnrdqAQfHMGd4Z/fa+oRUnRhBb /AEQ== X-Gm-Message-State: AIkVDXL0UfdEnMsP4ExihanBM6kIgIY0OzwlnZB2zyzTLwubMwCZET7L4zTaaW+ST5SwYw== X-Received: by 10.202.197.207 with SMTP id v198mr390626oif.71.1485921401875; Tue, 31 Jan 2017 19:56:41 -0800 (PST) Received: from Zany.attlocal.net (76-225-50-170.lightspeed.bkfdca.sbcglobal.net. [76.225.50.170]) by smtp.gmail.com with ESMTPSA id r131sm10036509oib.25.2017.01.31.19.56.40 for (version=TLS1_2 cipher=ECDHE-RSA-AES128-SHA bits=128/128); Tue, 31 Jan 2017 19:56:41 -0800 (PST) From: Thomas Turner To: ffmpeg-devel@ffmpeg.org Date: Tue, 31 Jan 2017 19:56:39 -0800 Message-Id: <1485921399-7462-1-git-send-email-thomastdt@googlemail.com> X-Mailer: git-send-email 1.9.1 Subject: [FFmpeg-devel] [PATCH] avutil/tests/lfg.c: added proper normality test X-BeenThere: ffmpeg-devel@ffmpeg.org X-Mailman-Version: 2.1.20 Precedence: list List-Id: FFmpeg development discussions and patches List-Unsubscribe: , List-Archive: List-Post: List-Help: List-Subscribe: , Reply-To: FFmpeg development discussions and patches MIME-Version: 1.0 Errors-To: ffmpeg-devel-bounces@ffmpeg.org Sender: "ffmpeg-devel" 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 --- 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