int cv::meanShift( InputArray _probImage, Rect& window, TermCriteria criteria ) { CvConnectedComp comp; Mat probImage = _probImage.getMat(); CvMat c_probImage = probImage; int iters = cvMeanShift(&c_probImage, window, (CvTermCriteria)criteria, &comp ); window = comp.rect; return iters; }
cvMeanShift( const void* imgProb, CvRect windowIn, CvTermCriteria criteria, CvConnectedComp* comp ) { CvMoments moments; int i = 0, eps; CvMat stub, *mat = (CvMat*)imgProb; CvMat cur_win; CvRect cur_rect = windowIn; if( comp ) comp->rect = windowIn; moments.m00 = moments.m10 = moments.m01 = 0; mat = cvGetMat( mat, &stub ); if( CV_MAT_CN( mat->type ) > 1 ) CV_Error( CV_BadNumChannels, cvUnsupportedFormat ); if( windowIn.height <= 0 || windowIn.width <= 0 ) CV_Error( CV_StsBadArg, "Input window has non-positive sizes" ); windowIn = cv::Rect(windowIn) & cv::Rect(0, 0, mat->cols, mat->rows); criteria = cvCheckTermCriteria( criteria, 1., 100 ); eps = cvRound( criteria.epsilon * criteria.epsilon ); for( i = 0; i < criteria.max_iter; i++ ) { int dx, dy, nx, ny; double inv_m00; cur_rect = cv::Rect(cur_rect) & cv::Rect(0, 0, mat->cols, mat->rows); if( cv::Rect(cur_rect) == cv::Rect() ) { cur_rect.x = mat->cols/2; cur_rect.y = mat->rows/2; } cur_rect.width = MAX(cur_rect.width, 1); cur_rect.height = MAX(cur_rect.height, 1); cvGetSubRect( mat, &cur_win, cur_rect ); cvMoments( &cur_win, &moments ); /* Calculating center of mass */ if( fabs(moments.m00) < DBL_EPSILON ) break; inv_m00 = moments.inv_sqrt_m00*moments.inv_sqrt_m00; dx = cvRound( moments.m10 * inv_m00 - windowIn.width*0.5 ); dy = cvRound( moments.m01 * inv_m00 - windowIn.height*0.5 ); nx = cur_rect.x + dx; ny = cur_rect.y + dy; if( nx < 0 ) nx = 0; else if( nx + cur_rect.width > mat->cols ) nx = mat->cols - cur_rect.width; if( ny < 0 ) ny = 0; else if( ny + cur_rect.height > mat->rows ) ny = mat->rows - cur_rect.height; dx = nx - cur_rect.x; dy = ny - cur_rect.y; cur_rect.x = nx; cur_rect.y = ny; /* Check for coverage centers mass & window */ if( dx*dx + dy*dy < eps ) break; } if( comp ) { comp->rect = cur_rect; comp->area = (float)moments.m00; } return i; }
typedef struct CvMoments { double m00, m10, m01, m20, m11, m02, m30, m21, m12, m03; /* spatial moments */ double mu20, mu11, mu02, mu30, mu21, mu12, mu03; /* central moments */ double inv_sqrt_m00; /* m00 != 0 ? 1/sqrt(m00) : 0 */ } CvMoments;
CV_IMPL void cvMoments( const void* array, CvMoments* moments, int binary ) { const int TILE_SIZE = 32; int type, depth, cn, coi = 0; CvMat stub, *mat = (CvMat*)array; CvMomentsInTileFunc func = 0; CvContour contourHeader; CvSeq* contour = 0; CvSeqBlock block; double buf[TILE_SIZE*TILE_SIZE]; uchar nzbuf[TILE_SIZE*TILE_SIZE]; if( CV_IS_SEQ( array )) { contour = (CvSeq*)array; if( !CV_IS_SEQ_POINT_SET( contour )) CV_Error( CV_StsBadArg, "The passed sequence is not a valid contour" ); } if( !moments ) CV_Error( CV_StsNullPtr, "" ); memset( moments, 0, sizeof(*moments)); if( !contour ) { mat = cvGetMat( mat, &stub, &coi ); type = CV_MAT_TYPE( mat->type ); if( type == CV_32SC2 || type == CV_32FC2 ) { contour = cvPointSeqFromMat( CV_SEQ_KIND_CURVE | CV_SEQ_FLAG_CLOSED, mat, &contourHeader, &block ); } } if( contour ) { icvContourMoments( contour, moments ); return; } type = CV_MAT_TYPE( mat->type ); depth = CV_MAT_DEPTH( type ); cn = CV_MAT_CN( type ); cv::Size size = cvGetMatSize( mat ); if( cn > 1 && coi == 0 ) CV_Error( CV_StsBadArg, "Invalid image type" ); if( size.width <= 0 || size.height <= 0 ) return; if( binary || depth == CV_8U ) func = momentsInTile<uchar, int, int>; else if( depth == CV_16U ) func = momentsInTile<ushort, int, int64>; else if( depth == CV_16S ) func = momentsInTile<short, int, int64>; else if( depth == CV_32F ) func = momentsInTile<float, double, double>; else if( depth == CV_64F ) func = momentsInTile<double, double, double>; else CV_Error( CV_StsUnsupportedFormat, "" ); cv::Mat src0(mat); for( int y = 0; y < size.height; y += TILE_SIZE ) { cv::Size tileSize; tileSize.height = std::min(TILE_SIZE, size.height - y); for( int x = 0; x < size.width; x += TILE_SIZE ) { tileSize.width = std::min(TILE_SIZE, size.width - x); cv::Mat src(src0, cv::Rect(x, y, tileSize.width, tileSize.height)); if( coi > 0 ) { cv::Mat tmp(tileSize, depth, buf); int pairs[] = {coi-1, 0}; cv::mixChannels(&src, 1, &tmp, 1, pairs, 1); src = tmp; } if( binary ) { cv::Mat tmp(tileSize, CV_8U, nzbuf); cv::compare( src, 0, tmp, CV_CMP_NE ); src = tmp; } double mom[10]; func( src, mom ); if(binary) { double s = 1./255; for( int k = 0; k < 10; k++ ) mom[k] *= s; } double xm = x * mom[0], ym = y * mom[0]; // accumulate moments computed in each tile // + m00 ( = m00' ) moments->m00 += mom[0]; // + m10 ( = m10' + x*m00' ) moments->m10 += mom[1] + xm; // + m01 ( = m01' + y*m00' ) moments->m01 += mom[2] + ym; // + m20 ( = m20' + 2*x*m10' + x*x*m00' ) moments->m20 += mom[3] + x * (mom[1] * 2 + xm); // + m11 ( = m11' + x*m01' + y*m10' + x*y*m00' ) moments->m11 += mom[4] + x * (mom[2] + ym) + y * mom[1]; // + m02 ( = m02' + 2*y*m01' + y*y*m00' ) moments->m02 += mom[5] + y * (mom[2] * 2 + ym); // + m30 ( = m30' + 3*x*m20' + 3*x*x*m10' + x*x*x*m00' ) moments->m30 += mom[6] + x * (3. * mom[3] + x * (3. * mom[1] + xm)); // + m21 ( = m21' + x*(2*m11' + 2*y*m10' + x*m01' + x*y*m00') + y*m20') moments->m21 += mom[7] + x * (2 * (mom[4] + y * mom[1]) + x * (mom[2] + ym)) + y * mom[3]; // + m12 ( = m12' + y*(2*m11' + 2*x*m01' + y*m10' + x*y*m00') + x*m02') moments->m12 += mom[8] + y * (2 * (mom[4] + x * mom[2]) + y * (mom[1] + xm)) + x * mom[5]; // + m03 ( = m03' + 3*y*m02' + 3*y*y*m01' + y*y*y*m00' ) moments->m03 += mom[9] + y * (3. * mom[5] + y * (3. * mom[2] + ym)); } } icvCompleteMomentState( moments ); }
参考:
https://blog.csdn.net/gdfsg/article/details/51016460