基于openCV的图像椭圆检测源代码（未使用库函数）

#include "stdafx.h" #include "stdio.h" #include "stdlib.h" #include "conio.h" #include "iostream" #include "math.h" #include "time.h" using namespace std; #include "cv.h" #include "highgui.h" #include "cxcore.h" #include "cvaux.h" using namespace cv; const int Bubblelength = 10000; //数组长度(防止数据溢出) const double pi = 3.141592653; //定义pi const double Threshold = 0.8; //投票阈值(判断此处是否存在椭圆) const double smallestAxis = 15; //检测椭圆大小阈值(只能检测半短轴大于smallestAxis像素的椭圆) int obtainMedian( Mat Input, int Size ); void obtain5parameters( CvMat *P, int Psize, double &a, double &b, double &c, double &d, double &e ); void votingFilter( CvMat *InMat, Mat ImageMat, int MatSize, CvSize ImageSize, Mat &bestOutput, double T, double sAxis ); void findRealEllipses( Mat BestM, int BestMsize, double Th, Mat &RealP, int &RealNum ); void main() { long begin,end; //估计算法运行时间 begin = clock(); //算法开始时间 int columnnum = 0; int contourlength,NumofrealEllipses; //轮廓长度，真实椭圆数量 char *fileName = "5.jpg"; //文件名 CvPoint center; //椭圆圆心坐标 CvSize axis; //长短半轴 CvSeq* contours; IplImage *originalImage = cvLoadImage( fileName, CV_LOAD_IMAGE_COLOR ); //将彩色图片导入 IplImage *ImageRGB = cvCreateImage(cvGetSize(originalImage), 8, 1 ); IplImage *SmoothyImage = cvCreateImage(cvGetSize(originalImage), 8, 1 ); IplImage *BImage = cvCreateImage(cvGetSize(originalImage), 8, 1 ); //创建8位1通道图像，用于存放二值图片 CvSize Imagesize = cvGetSize( originalImage ); Mat BImageMat = Mat::zeros( Imagesize.height, Imagesize.width, CV_64FC1 ); //存放二值化后图片像素矩阵 CvMemStorage *storage = cvCreateMemStorage(0); //创建内存存储器，用来存储返回的轮廓 cvConvertImage( originalImage, ImageRGB, 0 ); //将载入的图片转化为8位单通道图像放入ImageRGB cvSmooth( ImageRGB, SmoothyImage, CV_GAUSSIAN, 3, 3, 0 ); //平滑处理放入SmoothyImage cvCanny( SmoothyImage, BImage, 100, 255, 3); //取轮廓，获得二值图像 for ( int i = 0; i < Imagesize.height; i++ ) { for ( int j = 0; j < Imagesize.width; j ++) { BImageMat.at<double>(i,j) = cvGet2D(BImage,i,j).val[0]; } }//生成图像像素矩阵 cvNamedWindow( "Ellipse", 1 ); cvShowImage( "Ellipse", BImage ); //显示二值化后图片 cvWaitKey(20); int numberOfContours = cvFindContours( BImage, storage, &contours, sizeof(CvContour), CV_RETR_LIST, CV_CHAIN_APPROX_SIMPLE, cvPoint(0,0) );//获取轮廓 CvMat *fiveparmeters = cvCreateMat( 5, numberOfContours, CV_64FC1 ); //椭圆参数矩阵 for ( ; contours!=0; contours = contours -> h_next ) { //对每条轮廓逐个扫描输出参数向量 double A,B,C,D,E,p,q,semimajor,semiminor,theta,k; //定义椭圆参数及一些过度变量 contourlength = contours -> total; //获取每条轮廓的长度(每条轮廓所包含点的数量) if( contourlength >10) { //舍弃长度小于10的轮廓 CvMat *MatPinXY = cvCreateMat( 2, contourlength, CV_64FC1 ); for( int i = 0 ; i < contourlength; i++ ) { CvPoint *points = (CvPoint*)cvGetSeqElem(contours, i); cvmSet( MatPinXY, 0, i, points -> y); cvmSet( MatPinXY, 1, i, points -> x); }//将轮廓坐标写入MatPinXY(注意xy交换位置，转换为标准坐标系，方便计算) obtain5parameters( MatPinXY, contourlength, A, B, C, D, E ); //每条轮廓输出一个椭圆参数矩阵 } else { A = 0; B = 0; C = 0; D = 0; E = 0; } //对于不满足要求的轮廓，将椭圆参数向量定义为零 if( (B*B) < (4*A*C) ) { //参数向量满足椭圆约束判断 p = ((B*D-2*A*E)/(4*A*C-B*B)); q = ((B*E-2*C*D)/(4*A*C-B*B)); //求出圆心坐标(已经转换为图像坐标系) k = A*q*q+B*q*p+C*p*p-1; //过渡变量 if ( k*(A+C-sqrt(B*B+pow(A-C,2))) > 0 && p > 0 && q > 0) { //再次筛选 if( A < C) { //此时主轴长度大于副轴 semimajor = (sqrt((2*k)/(A+C-sqrt(B*B+pow((A-C),2))))); semiminor = (sqrt((2*k)/(A+C+sqrt(B*B+pow((A-C),2))))); } else { semimajor = (sqrt((2*k)/(A+C+sqrt(B*B+pow((A-C),2))))); semiminor = (sqrt((2*k)/(A+C-sqrt(B*B+pow((A-C),2))))); } if( (A-C)==0 ) { //偏角垂直判断 theta = 0.25*pi; } else { theta = 0.5*atan(B/(A-C)); } if( p < originalImage->width && q < originalImage->height) { //圆心坐标在图片内(再次筛选) cvmSet( fiveparmeters, 0, columnnum, p); cvmSet( fiveparmeters, 1, columnnum, q); cvmSet( fiveparmeters, 2, columnnum, semimajor); cvmSet( fiveparmeters, 3, columnnum, semiminor); cvmSet( fiveparmeters, 4, columnnum, theta); columnnum = columnnum + 1; }//输出符合要求的参数向量，并给出参数向量个数columnnum else {} } else{} } else{} } Mat best5Parameters = Mat::zeros( 5, columnnum, CV_64FC1 ); //定义最佳参数矩阵(这里列数取最值) Mat realEllipses = Mat::zeros( 5, columnnum, CV_64FC1 ); //定义真实椭圆参数矩阵(这里列数取最值) votingFilter( fiveparmeters, BImageMat, columnnum, Imagesize, best5Parameters, Threshold, smallestAxis );//投票滤波，输出最佳参数 findRealEllipses( best5Parameters, columnnum, smallestAxis, realEllipses, NumofrealEllipses ); end = clock(); //算法结束时间 cout << "此次计算耗时:"<< end-begin << "ms" << endl; //打印算法耗时 cout << endl; cout << "共检测到" << NumofrealEllipses << "个椭圆" << endl; cout << endl; for( int i =0; i < NumofrealEllipses; i++ ) { center.x = (int)realEllipses.at<double>(0,i); center.y = (int)realEllipses.at<double>(1,i); axis.width = (int)realEllipses.at<double>(2,i); axis.height = (int)realEllipses.at<double>(3,i); cvEllipse( originalImage, center, axis, 90-180*realEllipses.at<double>(4,i)/pi, 0, 360, CV_RGB(255, 0, 0), 1, CV_AA, 0); cout << "第" << i+1 << "个椭圆参数:" << endl; cout << "圆心坐标:(" << center.x << "," << center.y << ")" << endl; cout << "长短半轴:" << axis.width << " " << axis.height << endl; cout << "偏角" << 90-180*realEllipses.at<double>(4,i)/pi << "度" << endl; cout << endl; } cvNamedWindow( "Final", 1 ); cvShowImage( "Final", originalImage ); cvWaitKey(0); cvReleaseMemStorage( &storage ); cvReleaseMat( &fiveparmeters ); cvDestroyWindow( "Ellipse" ); cvDestroyWindow( "Final" ); } //获取每个轮廓对应的椭圆参数向量// void obtain5parameters( CvMat *P, int Psize, double &a, double &b, double &c, double &d, double &e ) { int counter = Psize/5; //循环计算次数(即随机取点的次数) int r,bestParaNum; //定义随机数r，最佳参数列向量标号(即在椭圆参数矩阵中的最佳参数向量) double ra,rb,rc,rd,re; //参数向量(并非最佳参数向量) Mat Para = Mat::zeros( 5, counter, CV_64FC1 ); //存放参数矩阵 Mat ParaD = Mat::zeros( 1, counter, CV_64FC1 ); //中值矩阵(存放每个参数向量对应的中值) Mat M = Mat::zeros( 5, 5, CV_64FC1 ); //线性方程组的左边矩阵 Mat X = Mat::zeros( 5, 1, CV_64FC1 ); //线性方程组的解向量 Mat Y = Mat::zeros( 5, 1, CV_64FC1 ); //线性方程组的右边向量 Mat MatD = Mat::zeros( 1, Psize, CV_64FC1 ); //存放每个点对应输出 for( int i = 0; i < 5; i++ ) { Y.at<double>(i,0) = -1; }//定义线性方程组右边矩阵 for( int i = 0; i < counter; i++ ) { for( int j = 0; j < 5; j++ ) { srand((i+1)*(j+1)); //随机数产生种子 r = rand()%Psize; //产生随机数(随机取点时使用) M.at<double>(j,0) = pow( cvmGet( P, 0, r ), 2 ); M.at<double>(j,1) = cvmGet( P, 1, r )*cvmGet( P, 0, r ); M.at<double>(j,2) = pow ( cvmGet( P, 1, r ), 2 ); M.at<double>(j,3) = cvmGet( P, 0, r ); M.at<double>(j,4) = cvmGet( P, 1, r ); }//通过随机采样的5个点求出线性方程组左边矩阵 Mat Mt = M.t(); Mat MtM = Mt*M; Mat Min = MtM.inv(); X = Min*Mt*Y;//求方程组的解向量(最小二乘法) for( int k = 0; k < 5; k++ ) { Para.at<double>(k,i) = X.at<double>(k,0); }//将参数向量放入参数矩阵 }//求出参�