在网上看了一篇关于cvHoughCircle函数的解析，转载过来看看。学习Opencv，不少人只能调用函数而不能作出适配或者修改，这是一种遗憾，貌似很少有这方面的博客，不知道有人推荐一些关于OpenCV源代码解析的资料没有？
一、基本原理
opencv中圆识别的基本原理如下：
1、canny算子求图像的单像素二值化边缘
2、假设我们需要找半径为R的所有圆，则对于边缘图中的每一个边缘点，该边缘点的切线的法线方向上(正负两个方向)，寻找到该边缘点距离为R的点，将该点的计数加1(初始化所有点的计数都是0)
3、找到计数值大于门限值的点，即圆心所在的点
二、代码分析
代码在/modules\imgproc\src\hough.cpp文件icvHoughCirclesGradient函数中
- static void
- icvHoughCirclesGradient( CvMat* img, float dp, float min_dist,
- int min_radius, int max_radius,
- int canny_threshold, int acc_threshold,
- CvSeq* circles, int circles_max )
- {
- //参数：
- //img: 输入图像
- //dp: 识别精度,1.0表示按照原图精度
- //min_dist: 圆心点位置识别精度
- //min_radius: 所需要找的圆的最小半径
- //max_radius：所需要找的圆的最大半径
- //canny_threshold：canny算子的高阀值
- //acc_threshold：累加器阀值，计数大于改阀值的点即被认为是可能的圆心
- //circles: 保存找到的符合条件的所有圆
- //circles_max: 最多需要的找到的圆的个数
- const int SHIFT = 10, ONE = 1 << SHIFT, R_THRESH = 30;
- cv::Ptr<CvMat> dx, dy;
- cv::Ptr<CvMat> edges, accum, dist_buf;
- std::vector<int> sort_buf;
- cv::Ptr<CvMemStorage> storage;
- int x, y, i, j, k, center_count, nz_count;
- float min_radius2 = (float)min_radius*min_radius;
- float max_radius2 = (float)max_radius*max_radius;
- int rows, cols, arows, acols;
- int astep, *adata;
- float* ddata;
- CvSeq *nz, *centers;
- float idp, dr;
- CvSeqReader reader;
- //canny算子求单像素二值化边缘，保存在edges变量中
- edges = cvCreateMat( img->rows, img->cols, CV_8UC1 );
- cvCanny( img, edges, MAX(canny_threshold/2,1), canny_threshold, 3 );
- //sobel算子求水平和垂直方向的边缘，用于计算边缘点的法线方向
- dx = cvCreateMat( img->rows, img->cols, CV_16SC1 );
- dy = cvCreateMat( img->rows, img->cols, CV_16SC1 );
- cvSobel( img, dx, 1, 0, 3 );
- cvSobel( img, dy, 0, 1, 3 );
- //dp表示识别精度
- if( dp < 1.f )
- dp = 1.f;
- idp = 1.f/dp;
- //accum用作累加器，包含图像中每一个点的计数。图像中每一个点都有一个计数，点的计数表示每一个canny边缘点法线方向上，
- //到该点距离为R的边缘点的个数，初始化为0
- accum = cvCreateMat( cvCeil(img->rows*idp)+2, cvCeil(img->cols*idp)+2, CV_32SC1 );
- cvZero(accum);
- storage = cvCreateMemStorage();
- nz = cvCreateSeq( CV_32SC2, sizeof(CvSeq), sizeof(CvPoint), storage );
- //centers用于保存可能的圆心点
- centers = cvCreateSeq( CV_32SC1, sizeof(CvSeq), sizeof(int), storage );
- rows = img->rows;
- cols = img->cols;
- arows = accum->rows - 2;
- acols = accum->cols - 2;
- adata = accum->data.i;
- astep = accum->step/sizeof(adata[0]);
- //以下这个循环用于获取所有可能的圆边缘点，存储在nz中，同时设置
- //累加器中的值
- for( y = 0; y < rows; y++ )
- {
- const uchar* edges_row = edges->data.ptr + y*edges->step;
- const short* dx_row = (const short*)(dx->data.ptr + y*dx->step);
- const short* dy_row = (const short*)(dy->data.ptr + y*dy->step);
- for( x = 0; x < cols; x++ )
- {
- float vx, vy;
- int sx, sy, x0, y0, x1, y1, r, k;
- CvPoint pt;
- vx = dx_row[x];
- vy = dy_row[x];
- if( !edges_row[x] || (vx == 0 && vy == 0) )
- continue;
- float mag = sqrt(vx*vx+vy*vy);
- assert( mag >= 1 );
- //sx表示cos, sy表示sin
- sx = cvRound((vx*idp)*ONE/mag);
- sy = cvRound((vy*idp)*ONE/mag);
- x0 = cvRound((x*idp)*ONE);
- y0 = cvRound((y*idp)*ONE);
- //循环两次表示需要计算两个方向，法线方向和法线的反方向
- for( k = 0; k < 2; k++ )
- {
- //半径方向的水平增量和垂直增量
- x1 = x0 + min_radius * sx;
- y1 = y0 + min_radius * sy;
- //在法线方向和反方向上，距离边缘点的距离为输入的最大半径和最小半径范围内找点
- //每找到一个点，该点的累加器计数就加1
- for( r = min_radius; r <= max_radius; x1 += sx, y1 += sy, r++ )
- {
- int x2 = x1 >> SHIFT, y2 = y1 >> SHIFT;
- if( (unsigned)x2 >= (unsigned)acols || (unsigned)y2 >= (unsigned)arows )
- break;
- adata[y2*astep + x2]++;
- }
- //反方向
- sx = -sx; sy = -sy;
- }
- //保存可能的圆边缘点
- pt.x = x; pt.y = y;
- cvSeqPush( nz, &pt );
- }
- }
- nz_count = nz->total;
- if( !nz_count )
- return;
- //累加器中，计数大于阀值的点，被认为可能的圆心点。因为计算各点计数过程中，距离有误差，所以
- //在与阀值比较时，该点计数先要与4邻域内的各个点的计数比较，最大者才能和阀值比较。可能的圆心
- //点保存在centers中
- for( y = 1; y < arows - 1; y++ )
- {
- for( x = 1; x < acols - 1; x++ )
- {
- int base = y*(acols+2) + x;
- if( adata[base] > acc_threshold &&
- adata[base] > adata[base-1] && adata[base] > adata[base+1] &&
- adata[base] > adata[base-acols-2] && adata[base] > adata[base+acols+2] )
- cvSeqPush(centers, &base);
- }
- }
- center_count = centers->total;
- if( !center_count )
- return;
- sort_buf.resize( MAX(center_count,nz_count) );
- //链表结构的certers转化成连续存储结构sort_buf
- cvCvtSeqToArray( centers, &sort_buf[0] );
- //经过icvHoughSortDescent32s函数后，以sort_buf中元素作为adata数组下标,
- //adata中的元素降序排列, 即adata[sort_buf[0]]是adata所有元素中最大的,
- //adata[sort_buf[center_count-1]]是所有元素中最小的
- icvHoughSortDescent32s( &sort_buf[0], center_count, adata );
- cvClearSeq( centers );
- //经过排序后的元素，重新以链表形式存储到centers中
- cvSeqPushMulti( centers, &sort_buf[0], center_count );
- dist_buf = cvCreateMat( 1, nz_count, CV_32FC1 );
- ddata = dist_buf->data.fl;
- dr = dp;
- min_dist = MAX( min_dist, dp );
- min_dist *= min_dist;
- //对于每一个可能的圆心点，计算所有边缘点到该圆心点的距离。由于centers中的
- //元素已经经过前面排序，所以累加器计数最大的可能圆心点最先进行下面的操作
- for( i = 0; i < centers->total; i++ )
- {
- int ofs = *(int*)cvGetSeqElem( centers, i );
- y = ofs/(acols+2) - 1;
- x = ofs - (y+1)*(acols+2) - 1;
- float cx = (float)(x*dp), cy = (float)(y*dp);
- float start_dist, dist_sum;
- float r_best = 0, c[3];
- int max_count = R_THRESH;
- //如果该可能的圆心点和已经确认的圆心点的距离小于阀值，则表示
- //这个圆心点和已经确认的圆心点是同一个点
- for( j = 0; j < circles->total; j++ )
- {
- float* c = (float*)cvGetSeqElem( circles, j );
- if( (c[0] - cx)*(c[0] - cx) + (c[1] - cy)*(c[1] - cy) < min_dist )
- break;
- }
- if( j < circles->total )
- continue;
- cvStartReadSeq( nz, &reader );
- //求所有边缘点到当前圆心点的距离，符合条件的距离值保存在ddata中
- for( j = k = 0; j < nz_count; j++ )
- {
- CvPoint pt;
- float _dx, _dy, _r2;
- CV_READ_SEQ_ELEM( pt, reader );
- _dx = cx - pt.x; _dy = cy - pt.y;
- _r2 = _dx*_dx + _dy*_dy;
- if(min_radius2 <= _r2 && _r2 <= max_radius2 )
- {
- ddata[k] = _r2;
- sort_buf[k] = k;
- k++;
- }
- }
- int nz_count1 = k, start_idx = nz_count1 - 1;
- if( nz_count1 == 0 )
- continue;
- dist_buf->cols = nz_count1;
- cvPow( dist_buf, dist_buf, 0.5 );
- //经过如下处理后，以sort_buf中元素作为ddata数组下标,ddata中的元素降序排列,
- //即ddata[sort_buf[0]]是ddata所有元素中最大的, ddata[sort_buf[nz_count1-1]]
- //是所有元素中最小的
- icvHoughSortDescent32s( &sort_buf[0], nz_count1, (int*)ddata );
- //对所有的距离值做处理，求出最可能圆半径值，max_count为到圆心的距离为最可能半径值的点的个数
- dist_sum = start_dist = ddata[sort_buf[nz_count1-1]];
- for( j = nz_count1 - 2; j >= 0; j-- )
- {
- float d = ddata[sort_buf[j]];
- if( d > max_radius )
- break;
- if( d - start_dist > dr )
- {
- float r_cur = ddata[sort_buf[(j + start_idx)/2]];
- if( (start_idx - j)*r_best >= max_count*r_cur ||
- (r_best < FLT_EPSILON && start_idx - j >= max_count) )
- {
- r_best = r_cur;
- max_count = start_idx - j;
- }
- start_dist = d;
- start_idx = j;
- dist_sum = 0;
- }
- dist_sum += d;
- }
- //max_count大于阀值，表示这几个边缘点构成一个圆
- if( max_count > R_THRESH )
- {
- c[0] = cx;
- c[1] = cy;
- c[2] = (float)r_best;
- cvSeqPush( circles, c );
- if( circles->total > circles_max )
- return;
- }
- }
- }
532
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
