1 // This file is part of OpenCV project.
2 // It is subject to the license terms in the LICENSE file found in the top-level directory
3 // of this distribution and at http://opencv.org/license.html.
4 #include "test_precomp.hpp"
6 namespace opencv_test { namespace {
8 class Core_RandTest : public cvtest::BaseTest
14 bool check_pdf(const Mat& hist, double scale, int dist_type,
15 double& refval, double& realval);
19 Core_RandTest::Core_RandTest()
23 static double chi2_p95(int n)
25 static float chi2_tab95[] = {
26 3.841f, 5.991f, 7.815f, 9.488f, 11.07f, 12.59f, 14.07f, 15.51f,
27 16.92f, 18.31f, 19.68f, 21.03f, 21.03f, 22.36f, 23.69f, 25.00f,
28 26.30f, 27.59f, 28.87f, 30.14f, 31.41f, 32.67f, 33.92f, 35.17f,
29 36.42f, 37.65f, 38.89f, 40.11f, 41.34f, 42.56f, 43.77f };
30 static const double xp = 1.64;
34 return chi2_tab95[n-1];
35 return n + sqrt((double)2*n)*xp + 0.6666666666666*(xp*xp - 1);
38 bool Core_RandTest::check_pdf(const Mat& hist, double scale,
39 int dist_type, double& refval, double& realval)
41 Mat hist0(hist.size(), CV_32F);
42 const int* H = hist.ptr<int>();
43 float* H0 = hist0.ptr<float>();
44 int i, hsz = hist.cols;
47 for( i = 0; i < hsz; i++ )
49 CV_Assert( fabs(1./sum - scale) < FLT_EPSILON );
51 if( dist_type == CV_RAND_UNI )
53 float scale0 = (float)(1./hsz);
54 for( i = 0; i < hsz; i++ )
59 double sum2 = 0, r = (hsz-1.)/2;
60 double alpha = 2*sqrt(2.)/r, beta = -alpha*r;
61 for( i = 0; i < hsz; i++ )
63 double x = i*alpha + beta;
64 H0[i] = (float)exp(-x*x);
68 for( i = 0; i < hsz; i++ )
69 H0[i] = (float)(H0[i]*sum2);
73 for( i = 0; i < hsz; i++ )
76 double b = H[i]*scale;
78 chi2 += (a - b)*(a - b)/(a + b);
82 double chi2_pval = chi2_p95(hsz - 1 - (dist_type == CV_RAND_NORMAL ? 2 : 0));
83 refval = chi2_pval*0.01;
84 return realval <= refval;
87 void Core_RandTest::run( int )
89 static int _ranges[][2] =
90 {{ 0, 256 }, { -128, 128 }, { 0, 65536 }, { -32768, 32768 },
91 { -1000000, 1000000 }, { -1000, 1000 }, { -1000, 1000 }};
93 const int MAX_SDIM = 10;
94 const int N = 2000000;
95 const int maxSlice = 1000;
96 const int MAX_HIST_SIZE = 1000;
99 RNG& rng = ts->get_rng();
100 RNG tested_rng = theRNG();
101 test_case_count = 200;
103 for( int idx = 0; idx < test_case_count; idx++ )
105 progress = update_progress( progress, idx, test_case_count, 0 );
106 ts->update_context( this, idx, false );
108 int depth = cvtest::randInt(rng) % (CV_64F+1);
109 int c, cn = (cvtest::randInt(rng) % 4) + 1;
110 int type = CV_MAKETYPE(depth, cn);
111 int dist_type = cvtest::randInt(rng) % (CV_RAND_NORMAL+1);
119 bool do_sphere_test = dist_type == CV_RAND_UNI;
123 arr[0].create(1, SZ, type);
124 arr[1].create(1, SZ, type);
125 bool fast_algo = dist_type == CV_RAND_UNI && depth < CV_32F;
127 for( c = 0; c < cn; c++ )
130 if( dist_type == CV_RAND_UNI )
132 a = (int)(cvtest::randInt(rng) % (_ranges[depth][1] -
133 _ranges[depth][0])) + _ranges[depth][0];
136 b = (int)(cvtest::randInt(rng) % (_ranges[depth][1] -
137 _ranges[depth][0])) + _ranges[depth][0];
139 while( abs(a-b) <= 1 );
143 unsigned r = (unsigned)(b - a);
144 fast_algo = fast_algo && r <= 256 && (r & (r-1)) == 0;
145 hsz = min((unsigned)(b - a), (unsigned)MAX_HIST_SIZE);
146 do_sphere_test = do_sphere_test && b - a >= 100;
150 int vrange = _ranges[depth][1] - _ranges[depth][0];
151 int meanrange = vrange/16;
152 int mindiv = MAX(vrange/20, 5);
153 int maxdiv = MIN(vrange/8, 10000);
155 a = cvtest::randInt(rng) % meanrange - meanrange/2 +
156 (_ranges[depth][0] + _ranges[depth][1])/2;
157 b = cvtest::randInt(rng) % (maxdiv - mindiv) + mindiv;
158 hsz = min((unsigned)b*9, (unsigned)MAX_HIST_SIZE);
162 hist[c].create(1, hsz, CV_32S);
165 cv::RNG saved_rng = tested_rng;
166 int maxk = fast_algo ? 0 : 1;
167 for( k = 0; k <= maxk; k++ )
169 tested_rng = saved_rng;
170 int sz = 0, dsz = 0, slice;
171 for( slice = 0; slice < maxSlice && sz < SZ; slice++, sz += dsz )
173 dsz = slice+1 < maxSlice ? (int)(cvtest::randInt(rng) % (SZ - sz) + 1) : SZ - sz;
174 Mat aslice = arr[k].colRange(sz, sz + dsz);
175 tested_rng.fill(aslice, dist_type, A, B);
176 //printf("%d - %d\n", sz, sz + dsz);
180 if( maxk >= 1 && cvtest::norm(arr[0], arr[1], NORM_INF) > eps)
182 ts->printf( cvtest::TS::LOG, "RNG output depends on the array lengths (some generated numbers get lost?)" );
183 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
187 for( c = 0; c < cn; c++ )
189 const uchar* data = arr[0].ptr();
190 int* H = hist[c].ptr<int>();
191 int HSZ = hist[c].cols;
192 double minVal = dist_type == CV_RAND_UNI ? A[c] : A[c] - B[c]*4;
193 double maxVal = dist_type == CV_RAND_UNI ? B[c] : A[c] + B[c]*4;
194 double scale = HSZ/(maxVal - minVal);
195 double delta = -minVal*scale;
197 hist[c] = Scalar::all(0);
199 for( i = c; i < SZ*cn; i += cn )
201 double val = depth == CV_8U ? ((const uchar*)data)[i] :
202 depth == CV_8S ? ((const schar*)data)[i] :
203 depth == CV_16U ? ((const ushort*)data)[i] :
204 depth == CV_16S ? ((const short*)data)[i] :
205 depth == CV_32S ? ((const int*)data)[i] :
206 depth == CV_32F ? ((const float*)data)[i] :
207 ((const double*)data)[i];
208 int ival = cvFloor(val*scale + delta);
209 if( (unsigned)ival < (unsigned)HSZ )
214 else if( dist_type == CV_RAND_UNI )
216 if( (minVal <= val && val < maxVal) || (depth >= CV_32F && val == maxVal) )
218 H[ival < 0 ? 0 : HSZ-1]++;
228 if( dist_type == CV_RAND_UNI && W[c] != SZ )
230 ts->printf( cvtest::TS::LOG, "Uniform RNG gave values out of the range [%g,%g) on channel %d/%d\n",
232 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
235 if( dist_type == CV_RAND_NORMAL && W[c] < SZ*.90)
237 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",
238 A[c], B[c], A[c], B[c], c, cn);
239 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
242 double refval = 0, realval = 0;
244 if( !check_pdf(hist[c], 1./W[c], dist_type, refval, realval) )
246 ts->printf( cvtest::TS::LOG, "RNG failed Chi-square test "
247 "(got %g vs probable maximum %g) on channel %d/%d\n",
248 realval, refval, c, cn);
249 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
254 // Monte-Carlo test. Compute volume of SDIM-dimensional sphere
255 // inscribed in [-1,1]^SDIM cube.
258 int SDIM = cvtest::randInt(rng) % (MAX_SDIM-1) + 2;
259 int N0 = (SZ*cn/SDIM), n = 0;
261 const uchar* data = arr[0].ptr();
262 double scale[4], delta[4];
263 for( c = 0; c < cn; c++ )
265 scale[c] = 2./(B[c] - A[c]);
266 delta[c] = -A[c]*scale[c] - 1;
269 for( i = k = c = 0; i <= SZ*cn - SDIM; i++, k++, c++ )
271 double val = depth == CV_8U ? ((const uchar*)data)[i] :
272 depth == CV_8S ? ((const schar*)data)[i] :
273 depth == CV_16U ? ((const ushort*)data)[i] :
274 depth == CV_16S ? ((const short*)data)[i] :
275 depth == CV_32S ? ((const int*)data)[i] :
276 depth == CV_32F ? ((const float*)data)[i] : ((const double*)data)[i];
277 c &= c < cn ? -1 : 0;
278 val = val*scale[c] + delta[c];
288 double V = ((double)n/N0)*(1 << SDIM);
290 // the theoretically computed volume
292 double V0 = sdim + 1;
293 for( sdim += 2; sdim <= SDIM; sdim += 2 )
296 if( fabs(V - V0) > 0.3*fabs(V0) )
298 ts->printf( cvtest::TS::LOG, "RNG failed %d-dim sphere volume test (got %g instead of %g)\n",
300 ts->printf( cvtest::TS::LOG, "depth = %d, N0 = %d\n", depth, N0);
301 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
308 TEST(Core_Rand, quality) { Core_RandTest test; test.safe_run(); }
311 class Core_RandRangeTest : public cvtest::BaseTest
314 Core_RandRangeTest() {}
315 ~Core_RandRangeTest() {}
319 Mat a(Size(1280, 720), CV_8U, Scalar(20));
320 Mat af(Size(1280, 720), CV_32F, Scalar(20));
321 theRNG().fill(a, RNG::UNIFORM, -DBL_MAX, DBL_MAX);
322 theRNG().fill(af, RNG::UNIFORM, -DBL_MAX, DBL_MAX);
323 int n0 = 0, n255 = 0, nx = 0;
324 int nfmin = 0, nfmax = 0, nfx = 0;
326 for( int i = 0; i < a.rows; i++ )
327 for( int j = 0; j < a.cols; j++ )
329 int v = a.at<uchar>(i,j);
330 double vf = af.at<float>(i,j);
332 else if( v == 255 ) n255++;
334 if( vf < FLT_MAX*-0.999f ) nfmin++;
335 else if( vf > FLT_MAX*0.999f ) nfmax++;
338 CV_Assert( n0 > nx*2 && n255 > nx*2 );
339 CV_Assert( nfmin > nfx*2 && nfmax > nfx*2 );
343 TEST(Core_Rand, range) { Core_RandRangeTest test; test.safe_run(); }
346 TEST(Core_RNG_MT19937, regression)
349 int actual[61] = {0, };
350 const size_t length = (sizeof(actual) / sizeof(actual[0]));
351 for (int i = 0; i < 10000; ++i )
353 actual[(unsigned)(rng.next() ^ i) % length]++;
356 int expected[length] = {
357 177, 158, 180, 177, 160, 179, 143, 162,
358 177, 144, 170, 174, 165, 168, 168, 156,
359 177, 157, 159, 169, 177, 182, 166, 154,
360 144, 180, 168, 152, 170, 187, 160, 145,
361 139, 164, 157, 179, 148, 183, 159, 160,
362 196, 184, 149, 142, 162, 148, 163, 152,
363 168, 173, 160, 181, 172, 181, 155, 153,
364 158, 171, 138, 150, 150 };
366 for (size_t i = 0; i < length; ++i)
368 ASSERT_EQ(expected[i], actual[i]);
373 TEST(Core_Rand, Regression_Stack_Corruption)
375 int bufsz = 128; //enough for 14 doubles
376 AutoBuffer<uchar> buffer(bufsz);
378 cv::Mat_<cv::Point2d> x(2, 3, (cv::Point2d*)(buffer.data()+offset));
379 offset += x.total()*x.elemSize();
380 double& param1 = *(double*)(buffer.data()+offset);
381 offset += sizeof(double);
382 double& param2 = *(double*)(buffer.data()+offset);
383 param1 = -9; param2 = 2;
385 cv::theRNG().fill(x, cv::RNG::NORMAL, param1, param2);
387 ASSERT_EQ(param1, -9);
388 ASSERT_EQ(param2, 2);
392 class RandRowFillParallelLoopBody : public cv::ParallelLoopBody
395 RandRowFillParallelLoopBody(Mat& dst) : dst_(dst) {}
396 ~RandRowFillParallelLoopBody() {}
397 void operator()(const cv::Range& r) const
399 cv::RNG rng = cv::theRNG(); // copy state
400 for (int y = r.start; y < r.end; y++)
402 cv::theRNG() = cv::RNG(rng.state + y); // seed is based on processed row
403 cv::randu(dst_.row(y), Scalar(-100), Scalar(100));
405 // theRNG() state is changed here (but state collision has low probability, so we don't check this)
411 TEST(Core_Rand, parallel_for_stable_results)
413 cv::RNG rng = cv::theRNG(); // save rng state
414 Mat dst1(1000, 100, CV_8SC1);
415 parallel_for_(cv::Range(0, dst1.rows), RandRowFillParallelLoopBody(dst1));
417 cv::theRNG() = rng; // restore rng state
418 Mat dst2(1000, 100, CV_8SC1);
419 parallel_for_(cv::Range(0, dst2.rows), RandRowFillParallelLoopBody(dst2));
421 ASSERT_EQ(0, countNonZero(dst1 != dst2));