in modules/core/src/rand.cpp [481:728]
void RNG::fill( InputOutputArray _mat, int disttype,
InputArray _param1arg, InputArray _param2arg, bool saturateRange )
{
Mat mat = _mat.getMat(), _param1 = _param1arg.getMat(), _param2 = _param2arg.getMat();
int depth = mat.depth(), cn = mat.channels();
AutoBuffer<double> _parambuf;
int j, k, fast_int_mode = 0, smallFlag = 1;
RandFunc func = 0;
RandnScaleFunc scaleFunc = 0;
CV_Assert(_param1.channels() == 1 && (_param1.rows == 1 || _param1.cols == 1) &&
(_param1.rows + _param1.cols - 1 == cn || _param1.rows + _param1.cols - 1 == 1 ||
(_param1.size() == Size(1, 4) && _param1.type() == CV_64F && cn <= 4)));
CV_Assert( _param2.channels() == 1 &&
(((_param2.rows == 1 || _param2.cols == 1) &&
(_param2.rows + _param2.cols - 1 == cn || _param2.rows + _param2.cols - 1 == 1 ||
(_param1.size() == Size(1, 4) && _param1.type() == CV_64F && cn <= 4))) ||
(_param2.rows == cn && _param2.cols == cn && disttype == NORMAL)));
Vec2i* ip = 0;
Vec2d* dp = 0;
Vec2f* fp = 0;
DivStruct* ds = 0;
uchar* mean = 0;
uchar* stddev = 0;
bool stdmtx = false;
int n1 = (int)_param1.total();
int n2 = (int)_param2.total();
if( disttype == UNIFORM )
{
_parambuf.allocate(cn*8 + n1 + n2);
double* parambuf = _parambuf;
double* p1 = _param1.ptr<double>();
double* p2 = _param2.ptr<double>();
if( !_param1.isContinuous() || _param1.type() != CV_64F || n1 != cn )
{
Mat tmp(_param1.size(), CV_64F, parambuf);
_param1.convertTo(tmp, CV_64F);
p1 = parambuf;
if( n1 < cn )
for( j = n1; j < cn; j++ )
p1[j] = p1[j-n1];
}
if( !_param2.isContinuous() || _param2.type() != CV_64F || n2 != cn )
{
Mat tmp(_param2.size(), CV_64F, parambuf + cn);
_param2.convertTo(tmp, CV_64F);
p2 = parambuf + cn;
if( n2 < cn )
for( j = n2; j < cn; j++ )
p2[j] = p2[j-n2];
}
if( depth <= CV_32S )
{
ip = (Vec2i*)(parambuf + cn*2);
for( j = 0, fast_int_mode = 1; j < cn; j++ )
{
double a = std::min(p1[j], p2[j]);
double b = std::max(p1[j], p2[j]);
if( saturateRange )
{
a = std::max(a, depth == CV_8U || depth == CV_16U ? 0. :
depth == CV_8S ? -128. : depth == CV_16S ? -32768. : (double)INT_MIN);
b = std::min(b, depth == CV_8U ? 256. : depth == CV_16U ? 65536. :
depth == CV_8S ? 128. : depth == CV_16S ? 32768. : (double)INT_MAX);
}
ip[j][1] = cvCeil(a);
int idiff = ip[j][0] = cvFloor(b) - ip[j][1] - 1;
double diff = b - a;
fast_int_mode &= diff <= 4294967296. && (idiff & (idiff+1)) == 0;
if( fast_int_mode )
smallFlag &= idiff <= 255;
else
{
if( diff > INT_MAX )
ip[j][0] = INT_MAX;
if( a < INT_MIN/2 )
ip[j][1] = INT_MIN/2;
}
}
if( !fast_int_mode )
{
ds = (DivStruct*)(ip + cn);
for( j = 0; j < cn; j++ )
{
ds[j].delta = ip[j][1];
unsigned d = ds[j].d = (unsigned)(ip[j][0]+1);
int l = 0;
while(((uint64)1 << l) < d)
l++;
ds[j].M = (unsigned)(((uint64)1 << 32)*(((uint64)1 << l) - d)/d) + 1;
ds[j].sh1 = std::min(l, 1);
ds[j].sh2 = std::max(l - 1, 0);
}
}
func = randTab[fast_int_mode][depth];
}
else
{
double scale = depth == CV_64F ?
5.4210108624275221700372640043497e-20 : // 2**-64
2.3283064365386962890625e-10; // 2**-32
double maxdiff = saturateRange ? (double)FLT_MAX : DBL_MAX;
// for each channel i compute such dparam[0][i] & dparam[1][i],
// so that a signed 32/64-bit integer X is transformed to
// the range [param1.val[i], param2.val[i]) using
// dparam[1][i]*X + dparam[0][i]
if( depth == CV_32F )
{
fp = (Vec2f*)(parambuf + cn*2);
for( j = 0; j < cn; j++ )
{
fp[j][0] = (float)(std::min(maxdiff, p2[j] - p1[j])*scale);
fp[j][1] = (float)((p2[j] + p1[j])*0.5);
}
}
else
{
dp = (Vec2d*)(parambuf + cn*2);
for( j = 0; j < cn; j++ )
{
dp[j][0] = std::min(DBL_MAX, p2[j] - p1[j])*scale;
dp[j][1] = ((p2[j] + p1[j])*0.5);
}
}
func = randTab[0][depth];
}
CV_Assert( func != 0 );
}
else if( disttype == CV_RAND_NORMAL )
{
_parambuf.allocate(MAX(n1, cn) + MAX(n2, cn));
double* parambuf = _parambuf;
int ptype = depth == CV_64F ? CV_64F : CV_32F;
int esz = (int)CV_ELEM_SIZE(ptype);
if( _param1.isContinuous() && _param1.type() == ptype )
mean = _param1.ptr();
else
{
Mat tmp(_param1.size(), ptype, parambuf);
_param1.convertTo(tmp, ptype);
mean = (uchar*)parambuf;
}
if( n1 < cn )
for( j = n1*esz; j < cn*esz; j++ )
mean[j] = mean[j - n1*esz];
if( _param2.isContinuous() && _param2.type() == ptype )
stddev = _param2.ptr();
else
{
Mat tmp(_param2.size(), ptype, parambuf + cn);
_param2.convertTo(tmp, ptype);
stddev = (uchar*)(parambuf + cn);
}
if( n1 < cn )
for( j = n1*esz; j < cn*esz; j++ )
stddev[j] = stddev[j - n1*esz];
stdmtx = _param2.rows == cn && _param2.cols == cn;
scaleFunc = randnScaleTab[depth];
CV_Assert( scaleFunc != 0 );
}
else
CV_Error( CV_StsBadArg, "Unknown distribution type" );
const Mat* arrays[] = {&mat, 0};
uchar* ptr;
NAryMatIterator it(arrays, &ptr);
int total = (int)it.size, blockSize = std::min((BLOCK_SIZE + cn - 1)/cn, total);
size_t esz = mat.elemSize();
AutoBuffer<double> buf;
uchar* param = 0;
float* nbuf = 0;
if( disttype == UNIFORM )
{
buf.allocate(blockSize*cn*4);
param = (uchar*)(double*)buf;
if( ip )
{
if( ds )
{
DivStruct* p = (DivStruct*)param;
for( j = 0; j < blockSize*cn; j += cn )
for( k = 0; k < cn; k++ )
p[j + k] = ds[k];
}
else
{
Vec2i* p = (Vec2i*)param;
for( j = 0; j < blockSize*cn; j += cn )
for( k = 0; k < cn; k++ )
p[j + k] = ip[k];
}
}
else if( fp )
{
Vec2f* p = (Vec2f*)param;
for( j = 0; j < blockSize*cn; j += cn )
for( k = 0; k < cn; k++ )
p[j + k] = fp[k];
}
else
{
Vec2d* p = (Vec2d*)param;
for( j = 0; j < blockSize*cn; j += cn )
for( k = 0; k < cn; k++ )
p[j + k] = dp[k];
}
}
else
{
buf.allocate((blockSize*cn+1)/2);
nbuf = (float*)(double*)buf;
}
for( size_t i = 0; i < it.nplanes; i++, ++it )
{
for( j = 0; j < total; j += blockSize )
{
int len = std::min(total - j, blockSize);
if( disttype == CV_RAND_UNI )
func( ptr, len*cn, &state, param, smallFlag != 0 );
else
{
randn_0_1_32f(nbuf, len*cn, &state);
scaleFunc(nbuf, ptr, len, cn, mean, stddev, stdmtx);
}
ptr += len*esz;
}
}
}