blob: 51408e960536b6ee608a2ab481960d305a4aee01 [file] [log] [blame]
jvanverth@google.com5a90ada2013-02-08 17:13:09 +00001
2/*
3 * Copyright 2013 Google Inc.
4 *
5 * Use of this source code is governed by a BSD-style license that can be
6 * found in the LICENSE file.
7 */
8
9#include "Test.h"
10#include "SkRandom.h"
11#include "SkTSort.h"
12
jvanverth@google.com897f4622013-02-08 18:28:47 +000013static bool anderson_darling_test(double p[32]) {
jvanverth@google.com5a90ada2013-02-08 17:13:09 +000014 // Min and max Anderson-Darling values allowable for k=32
15 const double kADMin32 = 0.202; // p-value of ~0.1
16 const double kADMax32 = 3.89; // p-value of ~0.99
17
18 // sort p values
19 SkTQSort<double>(p, p + 31);
20
21 // and compute Anderson-Darling statistic to ensure these are uniform
22 double s = 0.0;
23 for(int k = 0; k < 32; k++) {
24 double v = p[k]*(1.0 - p[31-k]);
25 if (v < 1.0e-30) {
26 v = 1.0e-30;
27 }
28 s += (2.0*(k+1)-1.0)*log(v);
29 }
30 double a2 = -32.0 - 0.03125*s;
31
32 return (kADMin32 < a2 && a2 < kADMax32);
33}
34
jvanverth@google.com897f4622013-02-08 18:28:47 +000035static bool chi_square_test(int bins[256], int e) {
jvanverth@google.com5a90ada2013-02-08 17:13:09 +000036 // Min and max chisquare values allowable
37 const double kChiSqMin256 = 206.3179; // probability of chance = 0.99 with k=256
38 const double kChiSqMax256 = 311.5603; // probability of chance = 0.01 with k=256
39
40 // compute chi-square
41 double chi2 = 0.0;
42 for (int j = 0; j < 256; ++j) {
43 double delta = bins[j] - e;
44 chi2 += delta*delta/e;
45 }
46
47 return (kChiSqMin256 < chi2 && chi2 < kChiSqMax256);
48}
49
50// Approximation to the normal distribution CDF
51// From Waissi and Rossin, 1996
jvanverth@google.com897f4622013-02-08 18:28:47 +000052static double normal_cdf(double z) {
jvanverth@google.com5a90ada2013-02-08 17:13:09 +000053 double t = ((-0.0004406*z*z* + 0.0418198)*z*z + 0.9)*z;
54 t *= -1.77245385091; // -sqrt(PI)
55 double p = 1.0/(1.0 + exp(t));
56
57 return p;
58}
59
jvanverth@google.com897f4622013-02-08 18:28:47 +000060static void test_random_byte(skiatest::Reporter* reporter, int shift) {
jvanverth@google.com5a90ada2013-02-08 17:13:09 +000061 int bins[256];
62 memset(bins, 0, sizeof(int)*256);
63
64 SkMWCRandom rand;
65 for (int i = 0; i < 256*10000; ++i) {
66 bins[(rand.nextU() >> shift) & 0xff]++;
67 }
skia.committer@gmail.com86267192013-02-09 07:05:02 +000068
jvanverth@google.com5a90ada2013-02-08 17:13:09 +000069 REPORTER_ASSERT(reporter, chi_square_test(bins, 10000));
70}
71
jvanverth@google.com897f4622013-02-08 18:28:47 +000072static void test_random_float(skiatest::Reporter* reporter) {
jvanverth@google.com5a90ada2013-02-08 17:13:09 +000073 int bins[256];
74 memset(bins, 0, sizeof(int)*256);
75
76 SkMWCRandom rand;
77 for (int i = 0; i < 256*10000; ++i) {
78 float f = rand.nextF();
79 REPORTER_ASSERT(reporter, 0.0f <= f && f < 1.0f);
80 bins[(int)(f*256.f)]++;
81 }
82 REPORTER_ASSERT(reporter, chi_square_test(bins, 10000));
83
84 double p[32];
85 for (int j = 0; j < 32; ++j) {
86 float f = rand.nextF();
87 REPORTER_ASSERT(reporter, 0.0f <= f && f < 1.0f);
88 p[j] = f;
89 }
90 REPORTER_ASSERT(reporter, anderson_darling_test(p));
91}
92
93// This is a test taken from tuftests by Marsaglia and Tsang. The idea here is that
94// we are using the random bit generated from a single shift position to generate
skia.committer@gmail.com86267192013-02-09 07:05:02 +000095// "strings" of 16 bits in length, shifting the string and adding a new bit with each
jvanverth@google.com5a90ada2013-02-08 17:13:09 +000096// iteration. We track the numbers generated. The ones that we don't generate will
skia.committer@gmail.com86267192013-02-09 07:05:02 +000097// have a normal distribution with mean ~24108 and standard deviation ~127. By
jvanverth@google.com5a90ada2013-02-08 17:13:09 +000098// creating a z-score (# of deviations from the mean) for one iteration of this step
99// we can determine its probability.
100//
skia.committer@gmail.com86267192013-02-09 07:05:02 +0000101// The original test used 26 bit strings, but is somewhat slow. This version uses 16
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000102// bits which is less rigorous but much faster to generate.
jvanverth@google.com897f4622013-02-08 18:28:47 +0000103static double test_single_gorilla(skiatest::Reporter* reporter, int shift) {
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000104 const int kWordWidth = 16;
105 const double kMean = 24108.0;
106 const double kStandardDeviation = 127.0;
107 const int kN = (1 << kWordWidth);
108 const int kNumEntries = kN >> 5; // dividing by 32
109 unsigned int entries[kNumEntries];
110
111 SkMWCRandom rand;
112 memset(entries, 0, sizeof(unsigned int)*kNumEntries);
113 // pre-seed our string value
114 int value = 0;
115 for (int i = 0; i < kWordWidth-1; ++i) {
116 value <<= 1;
117 unsigned int rnd = rand.nextU();
118 value |= ((rnd >> shift) & 0x1);
119 }
120
121 // now make some strings and track them
122 for (int i = 0; i < kN; ++i) {
123 value <<= 1;
124 unsigned int rnd = rand.nextU();
125 value |= ((rnd >> shift) & 0x1);
126
127 int index = value & (kNumEntries-1);
128 SkASSERT(index < kNumEntries);
jvanverth@google.com024e5232013-02-14 13:20:35 +0000129 int entry_shift = (value >> (kWordWidth-5)) & 0x1f;
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000130 entries[index] |= (0x1 << entry_shift);
131 }
132
133 // count entries
134 int total = 0;
135 for (int i = 0; i < kNumEntries; ++i) {
136 unsigned int entry = entries[i];
137 while (entry) {
138 total += (entry & 0x1);
139 entry >>= 1;
140 }
141 }
142
143 // convert counts to normal distribution z-score
skia.committer@gmail.com86267192013-02-09 07:05:02 +0000144 double z = ((kN-total)-kMean)/kStandardDeviation;
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000145
skia.committer@gmail.com86267192013-02-09 07:05:02 +0000146 // compute probability from normal distibution CDF
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000147 double p = normal_cdf(z);
148
jvanverth@google.com024e5232013-02-14 13:20:35 +0000149 REPORTER_ASSERT(reporter, 0.01 < p && p < 0.99);
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000150 return p;
151}
152
jvanverth@google.com897f4622013-02-08 18:28:47 +0000153static void test_gorilla(skiatest::Reporter* reporter) {
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000154
155 double p[32];
156 for (int bit_position = 0; bit_position < 32; ++bit_position) {
157 p[bit_position] = test_single_gorilla(reporter, bit_position);
158 }
skia.committer@gmail.com86267192013-02-09 07:05:02 +0000159
jvanverth@google.com024e5232013-02-14 13:20:35 +0000160 REPORTER_ASSERT(reporter, anderson_darling_test(p));
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000161}
162
jvanverth@google.com897f4622013-02-08 18:28:47 +0000163static void test_range(skiatest::Reporter* reporter) {
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000164 SkMWCRandom rand;
skia.committer@gmail.com86267192013-02-09 07:05:02 +0000165
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000166 // just to make sure we don't crash in this case
167 (void) rand.nextRangeU(0, 0xffffffff);
168
169 // check a case to see if it's uniform
170 int bins[256];
171 memset(bins, 0, sizeof(int)*256);
172 for (int i = 0; i < 256*10000; ++i) {
173 unsigned int u = rand.nextRangeU(17, 17+255);
174 REPORTER_ASSERT(reporter, 17 <= u && u <= 17+255);
175 bins[u - 17]++;
176 }
177
178 REPORTER_ASSERT(reporter, chi_square_test(bins, 10000));
179}
180
181static void TestRandom(skiatest::Reporter* reporter) {
182 // check uniform distributions of each byte in 32-bit word
183 test_random_byte(reporter, 0);
184 test_random_byte(reporter, 8);
185 test_random_byte(reporter, 16);
186 test_random_byte(reporter, 24);
187
188 test_random_float(reporter);
189
humper@google.com7cacbbd2013-02-08 18:44:27 +0000190 test_gorilla(reporter);
jvanverth@google.com5a90ada2013-02-08 17:13:09 +0000191
192 test_range(reporter);
193}
194
195#include "TestClassDef.h"
196DEFINE_TESTCLASS("Random", RandomTestClass, TestRandom)