上面一篇博客分析了HARRIS和ShiTomasi角点检测的源代码。而为了提取更准确的角点,OpenCV中提供了goodFeaturesToTrack()这个API函数,来获取更加准确的角点位置。这篇博客主要分析goodFeaturesToTrack()的源代码。
函数原型如下:
void cv::goodFeaturesToTrack( InputArray _image, OutputArray _corners, int maxCorners, double qualityLevel, double minDistance, InputArray _mask, int blockSize, int gradientSize, bool useHarrisDetector, double harrisK )
_image为输入的单通道图像;_corners为输出提取的角点坐标;maxCorners为设置的最大角点个数,程序中会按照角点强度降序排序,超过maxCorners的角点将被舍弃;qualityLevel为角点强度的阈值系数,如果检测出所有角点中的最大强度为max,则强度值<max*qualityLevel得角点会被舍弃;minDistance为角点与其邻域强角点之间的欧式距离,与邻域强角点距离小于minDistance的角点将被舍弃;_mask为设定的感兴趣区域,通常可不设置;blockSize是协方差矩阵滤波的窗口大小;gradientSize为sobel算子求微分的窗口的大小;useHarrisDetector为是否使用Harris角点检测;harrisK为Harris角点检测特征表达式中的常数k值。
接下来看源代码,首先根据useHarrisDetector这个flg值的设定选取不同的方法提取初始角点:
//提取初始角点
if( useHarrisDetector ) cornerHarris( image, eig, blockSize, gradientSize, harrisK ); else cornerMinEigenVal( image, eig, blockSize, gradientSize );
然后选取出所有角点特征值的最大值maxVal,并进行阈值化,将小于maxVal*qulityLevel的特征值舍弃掉(置为0),其余保持不变。并用3x3的方形膨胀核对特征值矩阵进行膨胀操作,膨胀的目的是使孤立的局部最大值扩大化,膨胀过后3x3邻域的非局部最大值点被局部最大值取代。膨胀过后的特征值图像为tmp,后面会用到。
//获取最大角点特征值 double maxVal = 0; minMaxLoc( eig, 0, &maxVal, 0, 0, _mask ); threshold( eig, eig, maxVal*qualityLevel, 0, THRESH_TOZERO ); dilate( eig, tmp, Mat());
取出局部最大值点,并将其地址保存到tmpCorners这个容器中。(val==tem_data[x]表示该点未被膨胀操作影响,就是局部极大值点)
// collect list of pointers to features - put them into temporary image Mat mask = _mask.getMat(); for( int y = 1; y < imgsize.height - 1; y++ ) { const float* eig_data = (const float*)eig.ptr(y); const float* tmp_data = (const float*)tmp.ptr(y); const uchar* mask_data = mask.data ? mask.ptr(y) : 0; for( int x = 1; x < imgsize.width - 1; x++ ) { float val = eig_data[x]; if( val != 0 && val == tmp_data[x] && (!mask_data || mask_data[x]) ) tmpCorners.push_back(eig_data + x); } }
然后将tmpCorners中指针指向的值进行降序排列,(其实是只改变指针指向的位置):
std::sort( tmpCorners.begin(), tmpCorners.end(), greaterThanPtr() );
接下来是不太好理解的地方,首先进行minDistance的判断。如果minDistance<1,也就是不进行minDistance的条件限制,直接保存前maxCorners个角点并返回;如果minDistance>=1,则执行if{}中的代码段。
首先弄清楚grid这个是什么东西,grid容器中元素个数是图像像素点个数的1/(cell_size*cell_size)倍,但是它是一个vector<vector<Point2f>>,也就是指grid中的每一个元素是一个vetcor<Point2f>(意思就是一个元素对应n个2维坐标点)。这个理解起来就是grid中相差一个像素,相当于原图像中相差minDistance个像素。原图像中以当前像素为中心的minDistance范围内的点在grid中其实视为一个像素的3x3邻域(一对多的关系)。
接下来的分析见注释:
if (minDistance >= 1) { // Partition the image into larger grids int w = image.cols; int h = image.rows; const int cell_size = cvRound(minDistance);
//grid中元素(vector)个数为grid_width*grid_height个,+cell_size-1的目的是保证能够覆盖所有的像素点 const int grid_width = (w + cell_size - 1) / cell_size; const int grid_height = (h + cell_size - 1) / cell_size; std::vector<std::vector<Point2f> > grid(grid_width*grid_height); minDistance *= minDistance; for( i = 0; i < total; i++ ) { int ofs = (int)((const uchar*)tmpCorners[i] - eig.ptr()); int y = (int)(ofs / eig.step); int x = (int)((ofs - y*eig.step)/sizeof(float)); //假设开始所有的角点是good bool good = true; //将原图中角点坐标归一化到mindistance邻域 int x_cell = x / cell_size; int y_cell = y / cell_size; //取归一化后点的邻域4个点坐标(不是4邻域,而是4个角),在grid中的坐标 int x1 = x_cell - 1; int y1 = y_cell - 1; int x2 = x_cell + 1; int y2 = y_cell + 1; // boundary check 边界判断 x1 = std::max(0, x1); y1 = std::max(0, y1); x2 = std::min(grid_width-1, x2); y2 = std::min(grid_height-1, y2); //遍历邻域4个点, grid中的邻域坐标 for( int yy = y1; yy <= y2; yy++ ) { for( int xx = x1; xx <= x2; xx++ ) {
//取出grid中一个元素对应的所有强角点坐标位置 std::vector <Point2f> &m = grid[yy*grid_width + xx]; //如果某元素对应的原图像角点容器中有已经保存的强角点,则需要进行距离判断。否则指定原图像中该角点就是强角点 if( m.size() ) {
//遍历其对应容器内的其他强角点,并依次判断原图像中当前角点与其邻域内其他强角点之间的欧式距离,如果欧式距离小于minDistance,则将当前角点标志置为good=false(抛弃),并跳出 for(j = 0; j < m.size(); j++) { float dx = x - m[j].x; float dy = y - m[j].y; if( dx*dx + dy*dy < minDistance ) { good = false; goto break_out; } } } } } break_out: //如果角点为good,则将该角点保存在当前grid中一个元素对应的容器中,同时保存在输出容器corners中,并累加计数器ncorners。由于已经进行过降序排序,前面保存的都是强角点。 if (good) { grid[y_cell*grid_width + x_cell].push_back(Point2f((float)x, (float)y)); corners.push_back(Point2f((float)x, (float)y)); ++ncorners; if( maxCorners > 0 && (int)ncorners == maxCorners ) break; } } } else { for( i = 0; i < total; i++ ) { int ofs = (int)((const uchar*)tmpCorners[i] - eig.ptr()); int y = (int)(ofs / eig.step); int x = (int)((ofs - y*eig.step)/sizeof(float)); corners.push_back(Point2f((float)x, (float)y)); ++ncorners; if( maxCorners > 0 && (int)ncorners == maxCorners ) break; } } Mat(corners).convertTo(_corners, _corners.fixedType() ? _corners.type() : CV_32F);
结束!