1 #include "test_precomp.hpp"
2
3 using namespace cv;
4 using namespace std;
5
6 class Core_RandTest : public cvtest::BaseTest
7 {
8 public:
9 Core_RandTest();
10 protected:
11 void run(int);
12 bool check_pdf(const Mat& hist, double scale, int dist_type,
13 double& refval, double& realval);
14 };
15
16
Core_RandTest()17 Core_RandTest::Core_RandTest()
18 {
19 }
20
chi2_p95(int n)21 static double chi2_p95(int n)
22 {
23 static float chi2_tab95[] = {
24 3.841f, 5.991f, 7.815f, 9.488f, 11.07f, 12.59f, 14.07f, 15.51f,
25 16.92f, 18.31f, 19.68f, 21.03f, 21.03f, 22.36f, 23.69f, 25.00f,
26 26.30f, 27.59f, 28.87f, 30.14f, 31.41f, 32.67f, 33.92f, 35.17f,
27 36.42f, 37.65f, 38.89f, 40.11f, 41.34f, 42.56f, 43.77f };
28 static const double xp = 1.64;
29 CV_Assert(n >= 1);
30
31 if( n <= 30 )
32 return chi2_tab95[n-1];
33 return n + sqrt((double)2*n)*xp + 0.6666666666666*(xp*xp - 1);
34 }
35
check_pdf(const Mat & hist,double scale,int dist_type,double & refval,double & realval)36 bool Core_RandTest::check_pdf(const Mat& hist, double scale,
37 int dist_type, double& refval, double& realval)
38 {
39 Mat hist0(hist.size(), CV_32F);
40 const int* H = hist.ptr<int>();
41 float* H0 = hist0.ptr<float>();
42 int i, hsz = hist.cols;
43
44 double sum = 0;
45 for( i = 0; i < hsz; i++ )
46 sum += H[i];
47 CV_Assert( fabs(1./sum - scale) < FLT_EPSILON );
48
49 if( dist_type == CV_RAND_UNI )
50 {
51 float scale0 = (float)(1./hsz);
52 for( i = 0; i < hsz; i++ )
53 H0[i] = scale0;
54 }
55 else
56 {
57 double sum2 = 0, r = (hsz-1.)/2;
58 double alpha = 2*sqrt(2.)/r, beta = -alpha*r;
59 for( i = 0; i < hsz; i++ )
60 {
61 double x = i*alpha + beta;
62 H0[i] = (float)exp(-x*x);
63 sum2 += H0[i];
64 }
65 sum2 = 1./sum2;
66 for( i = 0; i < hsz; i++ )
67 H0[i] = (float)(H0[i]*sum2);
68 }
69
70 double chi2 = 0;
71 for( i = 0; i < hsz; i++ )
72 {
73 double a = H0[i];
74 double b = H[i]*scale;
75 if( a > DBL_EPSILON )
76 chi2 += (a - b)*(a - b)/(a + b);
77 }
78 realval = chi2;
79
80 double chi2_pval = chi2_p95(hsz - 1 - (dist_type == CV_RAND_NORMAL ? 2 : 0));
81 refval = chi2_pval*0.01;
82 return realval <= refval;
83 }
84
run(int)85 void Core_RandTest::run( int )
86 {
87 static int _ranges[][2] =
88 {{ 0, 256 }, { -128, 128 }, { 0, 65536 }, { -32768, 32768 },
89 { -1000000, 1000000 }, { -1000, 1000 }, { -1000, 1000 }};
90
91 const int MAX_SDIM = 10;
92 const int N = 2000000;
93 const int maxSlice = 1000;
94 const int MAX_HIST_SIZE = 1000;
95 int progress = 0;
96
97 RNG& rng = ts->get_rng();
98 RNG tested_rng = theRNG();
99 test_case_count = 200;
100
101 for( int idx = 0; idx < test_case_count; idx++ )
102 {
103 progress = update_progress( progress, idx, test_case_count, 0 );
104 ts->update_context( this, idx, false );
105
106 int depth = cvtest::randInt(rng) % (CV_64F+1);
107 int c, cn = (cvtest::randInt(rng) % 4) + 1;
108 int type = CV_MAKETYPE(depth, cn);
109 int dist_type = cvtest::randInt(rng) % (CV_RAND_NORMAL+1);
110 int i, k, SZ = N/cn;
111 Scalar A, B;
112
113 double eps = 1.e-4;
114 if (depth == CV_64F)
115 eps = 1.e-7;
116
117 bool do_sphere_test = dist_type == CV_RAND_UNI;
118 Mat arr[2], hist[4];
119 int W[] = {0,0,0,0};
120
121 arr[0].create(1, SZ, type);
122 arr[1].create(1, SZ, type);
123 bool fast_algo = dist_type == CV_RAND_UNI && depth < CV_32F;
124
125 for( c = 0; c < cn; c++ )
126 {
127 int a, b, hsz;
128 if( dist_type == CV_RAND_UNI )
129 {
130 a = (int)(cvtest::randInt(rng) % (_ranges[depth][1] -
131 _ranges[depth][0])) + _ranges[depth][0];
132 do
133 {
134 b = (int)(cvtest::randInt(rng) % (_ranges[depth][1] -
135 _ranges[depth][0])) + _ranges[depth][0];
136 }
137 while( abs(a-b) <= 1 );
138 if( a > b )
139 std::swap(a, b);
140
141 unsigned r = (unsigned)(b - a);
142 fast_algo = fast_algo && r <= 256 && (r & (r-1)) == 0;
143 hsz = min((unsigned)(b - a), (unsigned)MAX_HIST_SIZE);
144 do_sphere_test = do_sphere_test && b - a >= 100;
145 }
146 else
147 {
148 int vrange = _ranges[depth][1] - _ranges[depth][0];
149 int meanrange = vrange/16;
150 int mindiv = MAX(vrange/20, 5);
151 int maxdiv = MIN(vrange/8, 10000);
152
153 a = cvtest::randInt(rng) % meanrange - meanrange/2 +
154 (_ranges[depth][0] + _ranges[depth][1])/2;
155 b = cvtest::randInt(rng) % (maxdiv - mindiv) + mindiv;
156 hsz = min((unsigned)b*9, (unsigned)MAX_HIST_SIZE);
157 }
158 A[c] = a;
159 B[c] = b;
160 hist[c].create(1, hsz, CV_32S);
161 }
162
163 cv::RNG saved_rng = tested_rng;
164 int maxk = fast_algo ? 0 : 1;
165 for( k = 0; k <= maxk; k++ )
166 {
167 tested_rng = saved_rng;
168 int sz = 0, dsz = 0, slice;
169 for( slice = 0; slice < maxSlice; slice++, sz += dsz )
170 {
171 dsz = slice+1 < maxSlice ? (int)(cvtest::randInt(rng) % (SZ - sz + 1)) : SZ - sz;
172 Mat aslice = arr[k].colRange(sz, sz + dsz);
173 tested_rng.fill(aslice, dist_type, A, B);
174 }
175 }
176
177 if( maxk >= 1 && cvtest::norm(arr[0], arr[1], NORM_INF) > eps)
178 {
179 ts->printf( cvtest::TS::LOG, "RNG output depends on the array lengths (some generated numbers get lost?)" );
180 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
181 return;
182 }
183
184 for( c = 0; c < cn; c++ )
185 {
186 const uchar* data = arr[0].ptr();
187 int* H = hist[c].ptr<int>();
188 int HSZ = hist[c].cols;
189 double minVal = dist_type == CV_RAND_UNI ? A[c] : A[c] - B[c]*4;
190 double maxVal = dist_type == CV_RAND_UNI ? B[c] : A[c] + B[c]*4;
191 double scale = HSZ/(maxVal - minVal);
192 double delta = -minVal*scale;
193
194 hist[c] = Scalar::all(0);
195
196 for( i = c; i < SZ*cn; i += cn )
197 {
198 double val = depth == CV_8U ? ((const uchar*)data)[i] :
199 depth == CV_8S ? ((const schar*)data)[i] :
200 depth == CV_16U ? ((const ushort*)data)[i] :
201 depth == CV_16S ? ((const short*)data)[i] :
202 depth == CV_32S ? ((const int*)data)[i] :
203 depth == CV_32F ? ((const float*)data)[i] :
204 ((const double*)data)[i];
205 int ival = cvFloor(val*scale + delta);
206 if( (unsigned)ival < (unsigned)HSZ )
207 {
208 H[ival]++;
209 W[c]++;
210 }
211 else if( dist_type == CV_RAND_UNI )
212 {
213 if( (minVal <= val && val < maxVal) || (depth >= CV_32F && val == maxVal) )
214 {
215 H[ival < 0 ? 0 : HSZ-1]++;
216 W[c]++;
217 }
218 else
219 {
220 putchar('^');
221 }
222 }
223 }
224
225 if( dist_type == CV_RAND_UNI && W[c] != SZ )
226 {
227 ts->printf( cvtest::TS::LOG, "Uniform RNG gave values out of the range [%g,%g) on channel %d/%d\n",
228 A[c], B[c], c, cn);
229 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
230 return;
231 }
232 if( dist_type == CV_RAND_NORMAL && W[c] < SZ*.90)
233 {
234 ts->printf( cvtest::TS::LOG, "Normal RNG gave too many values out of the range (%g+4*%g,%g+4*%g) on channel %d/%d\n",
235 A[c], B[c], A[c], B[c], c, cn);
236 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
237 return;
238 }
239 double refval = 0, realval = 0;
240
241 if( !check_pdf(hist[c], 1./W[c], dist_type, refval, realval) )
242 {
243 ts->printf( cvtest::TS::LOG, "RNG failed Chi-square test "
244 "(got %g vs probable maximum %g) on channel %d/%d\n",
245 realval, refval, c, cn);
246 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
247 return;
248 }
249 }
250
251 // Monte-Carlo test. Compute volume of SDIM-dimensional sphere
252 // inscribed in [-1,1]^SDIM cube.
253 if( do_sphere_test )
254 {
255 int SDIM = cvtest::randInt(rng) % (MAX_SDIM-1) + 2;
256 int N0 = (SZ*cn/SDIM), n = 0;
257 double r2 = 0;
258 const uchar* data = arr[0].ptr();
259 double scale[4], delta[4];
260 for( c = 0; c < cn; c++ )
261 {
262 scale[c] = 2./(B[c] - A[c]);
263 delta[c] = -A[c]*scale[c] - 1;
264 }
265
266 for( i = k = c = 0; i <= SZ*cn - SDIM; i++, k++, c++ )
267 {
268 double val = depth == CV_8U ? ((const uchar*)data)[i] :
269 depth == CV_8S ? ((const schar*)data)[i] :
270 depth == CV_16U ? ((const ushort*)data)[i] :
271 depth == CV_16S ? ((const short*)data)[i] :
272 depth == CV_32S ? ((const int*)data)[i] :
273 depth == CV_32F ? ((const float*)data)[i] : ((const double*)data)[i];
274 c &= c < cn ? -1 : 0;
275 val = val*scale[c] + delta[c];
276 r2 += val*val;
277 if( k == SDIM-1 )
278 {
279 n += r2 <= 1;
280 r2 = 0;
281 k = -1;
282 }
283 }
284
285 double V = ((double)n/N0)*(1 << SDIM);
286
287 // the theoretically computed volume
288 int sdim = SDIM % 2;
289 double V0 = sdim + 1;
290 for( sdim += 2; sdim <= SDIM; sdim += 2 )
291 V0 *= 2*CV_PI/sdim;
292
293 if( fabs(V - V0) > 0.3*fabs(V0) )
294 {
295 ts->printf( cvtest::TS::LOG, "RNG failed %d-dim sphere volume test (got %g instead of %g)\n",
296 SDIM, V, V0);
297 ts->printf( cvtest::TS::LOG, "depth = %d, N0 = %d\n", depth, N0);
298 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
299 return;
300 }
301 }
302 }
303 }
304
TEST(Core_Rand,quality)305 TEST(Core_Rand, quality) { Core_RandTest test; test.safe_run(); }
306
307
308 class Core_RandRangeTest : public cvtest::BaseTest
309 {
310 public:
Core_RandRangeTest()311 Core_RandRangeTest() {}
~Core_RandRangeTest()312 ~Core_RandRangeTest() {}
313 protected:
run(int)314 void run(int)
315 {
316 Mat a(Size(1280, 720), CV_8U, Scalar(20));
317 Mat af(Size(1280, 720), CV_32F, Scalar(20));
318 theRNG().fill(a, RNG::UNIFORM, -DBL_MAX, DBL_MAX);
319 theRNG().fill(af, RNG::UNIFORM, -DBL_MAX, DBL_MAX);
320 int n0 = 0, n255 = 0, nx = 0;
321 int nfmin = 0, nfmax = 0, nfx = 0;
322
323 for( int i = 0; i < a.rows; i++ )
324 for( int j = 0; j < a.cols; j++ )
325 {
326 int v = a.at<uchar>(i,j);
327 double vf = af.at<float>(i,j);
328 if( v == 0 ) n0++;
329 else if( v == 255 ) n255++;
330 else nx++;
331 if( vf < FLT_MAX*-0.999f ) nfmin++;
332 else if( vf > FLT_MAX*0.999f ) nfmax++;
333 else nfx++;
334 }
335 CV_Assert( n0 > nx*2 && n255 > nx*2 );
336 CV_Assert( nfmin > nfx*2 && nfmax > nfx*2 );
337 }
338 };
339
TEST(Core_Rand,range)340 TEST(Core_Rand, range) { Core_RandRangeTest test; test.safe_run(); }
341
342
TEST(Core_RNG_MT19937,regression)343 TEST(Core_RNG_MT19937, regression)
344 {
345 cv::RNG_MT19937 rng;
346 int actual[61] = {0, };
347 const size_t length = (sizeof(actual) / sizeof(actual[0]));
348 for (int i = 0; i < 10000; ++i )
349 {
350 actual[(unsigned)(rng.next() ^ i) % length]++;
351 }
352
353 int expected[length] = {
354 177, 158, 180, 177, 160, 179, 143, 162,
355 177, 144, 170, 174, 165, 168, 168, 156,
356 177, 157, 159, 169, 177, 182, 166, 154,
357 144, 180, 168, 152, 170, 187, 160, 145,
358 139, 164, 157, 179, 148, 183, 159, 160,
359 196, 184, 149, 142, 162, 148, 163, 152,
360 168, 173, 160, 181, 172, 181, 155, 153,
361 158, 171, 138, 150, 150 };
362
363 for (size_t i = 0; i < length; ++i)
364 {
365 ASSERT_EQ(expected[i], actual[i]);
366 }
367 }
368