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; 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);
179 if( maxk >= 1 && cvtest::norm(arr[0], arr[1], NORM_INF) > eps)
181 ts->printf( cvtest::TS::LOG, "RNG output depends on the array lengths (some generated numbers get lost?)" );
182 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
186 for( c = 0; c < cn; c++ )
188 const uchar* data = arr[0].ptr();
189 int* H = hist[c].ptr<int>();
190 int HSZ = hist[c].cols;
191 double minVal = dist_type == CV_RAND_UNI ? A[c] : A[c] - B[c]*4;
192 double maxVal = dist_type == CV_RAND_UNI ? B[c] : A[c] + B[c]*4;
193 double scale = HSZ/(maxVal - minVal);
194 double delta = -minVal*scale;
196 hist[c] = Scalar::all(0);
198 for( i = c; i < SZ*cn; i += cn )
200 double val = depth == CV_8U ? ((const uchar*)data)[i] :
201 depth == CV_8S ? ((const schar*)data)[i] :
202 depth == CV_16U ? ((const ushort*)data)[i] :
203 depth == CV_16S ? ((const short*)data)[i] :
204 depth == CV_32S ? ((const int*)data)[i] :
205 depth == CV_32F ? ((const float*)data)[i] :
206 ((const double*)data)[i];
207 int ival = cvFloor(val*scale + delta);
208 if( (unsigned)ival < (unsigned)HSZ )
213 else if( dist_type == CV_RAND_UNI )
215 if( (minVal <= val && val < maxVal) || (depth >= CV_32F && val == maxVal) )
217 H[ival < 0 ? 0 : HSZ-1]++;
227 if( dist_type == CV_RAND_UNI && W[c] != SZ )
229 ts->printf( cvtest::TS::LOG, "Uniform RNG gave values out of the range [%g,%g) on channel %d/%d\n",
231 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
234 if( dist_type == CV_RAND_NORMAL && W[c] < SZ*.90)
236 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",
237 A[c], B[c], A[c], B[c], c, cn);
238 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
241 double refval = 0, realval = 0;
243 if( !check_pdf(hist[c], 1./W[c], dist_type, refval, realval) )
245 ts->printf( cvtest::TS::LOG, "RNG failed Chi-square test "
246 "(got %g vs probable maximum %g) on channel %d/%d\n",
247 realval, refval, c, cn);
248 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
253 // Monte-Carlo test. Compute volume of SDIM-dimensional sphere
254 // inscribed in [-1,1]^SDIM cube.
257 int SDIM = cvtest::randInt(rng) % (MAX_SDIM-1) + 2;
258 int N0 = (SZ*cn/SDIM), n = 0;
260 const uchar* data = arr[0].ptr();
261 double scale[4], delta[4];
262 for( c = 0; c < cn; c++ )
264 scale[c] = 2./(B[c] - A[c]);
265 delta[c] = -A[c]*scale[c] - 1;
268 for( i = k = c = 0; i <= SZ*cn - SDIM; i++, k++, c++ )
270 double val = depth == CV_8U ? ((const uchar*)data)[i] :
271 depth == CV_8S ? ((const schar*)data)[i] :
272 depth == CV_16U ? ((const ushort*)data)[i] :
273 depth == CV_16S ? ((const short*)data)[i] :
274 depth == CV_32S ? ((const int*)data)[i] :
275 depth == CV_32F ? ((const float*)data)[i] : ((const double*)data)[i];
276 c &= c < cn ? -1 : 0;
277 val = val*scale[c] + delta[c];
287 double V = ((double)n/N0)*(1 << SDIM);
289 // the theoretically computed volume
291 double V0 = sdim + 1;
292 for( sdim += 2; sdim <= SDIM; sdim += 2 )
295 if( fabs(V - V0) > 0.3*fabs(V0) )
297 ts->printf( cvtest::TS::LOG, "RNG failed %d-dim sphere volume test (got %g instead of %g)\n",
299 ts->printf( cvtest::TS::LOG, "depth = %d, N0 = %d\n", depth, N0);
300 ts->set_failed_test_info( cvtest::TS::FAIL_INVALID_OUTPUT );
307 TEST(Core_Rand, quality) { Core_RandTest test; test.safe_run(); }
310 class Core_RandRangeTest : public cvtest::BaseTest
313 Core_RandRangeTest() {}
314 ~Core_RandRangeTest() {}
318 Mat a(Size(1280, 720), CV_8U, Scalar(20));
319 Mat af(Size(1280, 720), CV_32F, Scalar(20));
320 theRNG().fill(a, RNG::UNIFORM, -DBL_MAX, DBL_MAX);
321 theRNG().fill(af, RNG::UNIFORM, -DBL_MAX, DBL_MAX);
322 int n0 = 0, n255 = 0, nx = 0;
323 int nfmin = 0, nfmax = 0, nfx = 0;
325 for( int i = 0; i < a.rows; i++ )
326 for( int j = 0; j < a.cols; j++ )
328 int v = a.at<uchar>(i,j);
329 double vf = af.at<float>(i,j);
331 else if( v == 255 ) n255++;
333 if( vf < FLT_MAX*-0.999f ) nfmin++;
334 else if( vf > FLT_MAX*0.999f ) nfmax++;
337 CV_Assert( n0 > nx*2 && n255 > nx*2 );
338 CV_Assert( nfmin > nfx*2 && nfmax > nfx*2 );
342 TEST(Core_Rand, range) { Core_RandRangeTest test; test.safe_run(); }
345 TEST(Core_RNG_MT19937, regression)
348 int actual[61] = {0, };
349 const size_t length = (sizeof(actual) / sizeof(actual[0]));
350 for (int i = 0; i < 10000; ++i )
352 actual[(unsigned)(rng.next() ^ i) % length]++;
355 int expected[length] = {
356 177, 158, 180, 177, 160, 179, 143, 162,
357 177, 144, 170, 174, 165, 168, 168, 156,
358 177, 157, 159, 169, 177, 182, 166, 154,
359 144, 180, 168, 152, 170, 187, 160, 145,
360 139, 164, 157, 179, 148, 183, 159, 160,
361 196, 184, 149, 142, 162, 148, 163, 152,
362 168, 173, 160, 181, 172, 181, 155, 153,
363 158, 171, 138, 150, 150 };
365 for (size_t i = 0; i < length; ++i)
367 ASSERT_EQ(expected[i], actual[i]);
372 TEST(Core_Rand, Regression_Stack_Corruption)
374 int bufsz = 128; //enough for 14 doubles
375 AutoBuffer<uchar> buffer(bufsz);
377 cv::Mat_<cv::Point2d> x(2, 3, (cv::Point2d*)(buffer+offset)); offset += x.total()*x.elemSize();
378 double& param1 = *(double*)(buffer+offset); offset += sizeof(double);
379 double& param2 = *(double*)(buffer+offset); offset += sizeof(double);
380 param1 = -9; param2 = 2;
382 cv::theRNG().fill(x, cv::RNG::NORMAL, param1, param2);
384 ASSERT_EQ(param1, -9);
385 ASSERT_EQ(param2, 2);
389 class RandRowFillParallelLoopBody : public cv::ParallelLoopBody
392 RandRowFillParallelLoopBody(Mat& dst) : dst_(dst) {}
393 ~RandRowFillParallelLoopBody() {}
394 void operator()(const cv::Range& r) const
396 cv::RNG rng = cv::theRNG(); // copy state
397 for (int y = r.start; y < r.end; y++)
399 cv::theRNG() = cv::RNG(rng.state + y); // seed is based on processed row
400 cv::randu(dst_.row(y), Scalar(-100), Scalar(100));
402 // theRNG() state is changed here (but state collision has low probability, so we don't check this)
408 TEST(Core_Rand, parallel_for_stable_results)
410 cv::RNG rng = cv::theRNG(); // save rng state
411 Mat dst1(1000, 100, CV_8SC1);
412 parallel_for_(cv::Range(0, dst1.rows), RandRowFillParallelLoopBody(dst1));
414 cv::theRNG() = rng; // restore rng state
415 Mat dst2(1000, 100, CV_8SC1);
416 parallel_for_(cv::Range(0, dst2.rows), RandRowFillParallelLoopBody(dst2));
418 ASSERT_EQ(0, countNonZero(dst1 != dst2));