霍夫变换检测椭圆

2022-07-04 智能车 0 评论

霍夫变换（一）

//定义椭圆结构体，成员包含累积计数，中心点，长短半轴，倾斜角度 
typedef struct Ellipse {
    int accumulator;
    double xc, yc;
    double a, b;
    double angle;
} Ellipse;

//创建格子，传累加器列表和精度 
static std::vector<double> _create_bins(const std::vector<double>& acc, double binSize) {
	//把累加器列表中最大值取出，即计算得到的短半轴平方的所有值中最大的值 
    std::vector<double>::const_iterator maxPos = std::max_element(acc.begin(), acc.end());
    double maxVal = *maxPos;
    //要返回的格子列表 
    std::vector<double> v;
    //格子是从0开始，到最大值，间隔为精度 
    for (double d = 0; d < maxVal + binSize; d += binSize) {
        v.push_back(d);
    }
    return v;
}

//然后是计算直方图数据的函数，传入累加器列表和格子列表，返回就是累加器列表中的值落入相应的格子的计算的列表 
static std::vector<int> _histogram(const std::vector<double>& acc, const std::vector<double>& bins) {
	//直方图的数据列表，大小比格子的边界数少1 
    std::vector<int> hist(bins.size()-1);
    //直方图数据初始化为0 
    for (int i = 0; i < hist.size(); i++) {
        hist[i] = 0;
    }
    //遍历累加器列表中的数据，看落入了哪个格子 
    for (int i = 0; i < acc.size(); i++) {
        double a = acc[i];
        //看落入哪个格子 
        for (int j = 0; j < bins.size()-1; j++) {
        	//格子的左边值和右边值，除了最后一个格子，都是左闭右开 
            double left = bins[j];
            double right = bins[j+1];
            //如果值落入相应的格子，相应的直方图计算加1 
            if (a >= left && a < right) {
                hist[j]++;
                break;
            }
        }
        //如果为最后一个格子的边界，最后一个格子的计算也要加1，最后一个格子是左闭右也闭 
        if (a == bins[bins.size()-1]) {
            hist[bins.size()-2]++;
        }
    }
    return hist;
}

//定义算法的实现函数，返回检测到的椭圆列表
//参数：传入的边缘图，
//累计计数的阈值（即达到什么值才认为有效），
//精度（主要用来定义直方图的格子），
//轴的最小值，
//轴的最大值（如果为-1，则取图像宽高中小的那个一半），
//第一个轴的最大值（如果为-1，则无限制） 
static std::vector<Ellipse> _hough_ellipse(const cv::Mat& edges, int threshold, double accuracy,
 int minAxis, int maxAxis=-1, int maxFirst=-1) {
    //定义返回的列表 
	std::vector<Ellipse> results;
    
    //拿到图像的宽高 
    int numRows = edges.rows;
    int numCols = edges.cols;
    //边缘图像的非零点位置 
    std::vector<int> nonzeroXIndices, nonzeroYIndices;
    //把边缘图像的非零点找出来 
    for (int i = 0; i < numRows; i++) {
        for (int j = 0; j < numCols; j++) {
            if (edges.ptr<uchar>(i)[j] > 0) {
                nonzeroXIndices.push_back(j);
                nonzeroYIndices.push_back(i);
            }
        }
    }
    //总共有多少非零点 
    int numPixels = nonzeroXIndices.size();
    //累加器列表 
    std::vector<double> accumulator;
    //精度，即直方图的格子宽度，这里取平方，因为下面的算法计算短半轴时为了减少计算量也是用平方（不开方） 
    double binSize = accuracy*accuracy;
    //最大的半轴值，也取平方以避免开方 
    double maxAxisSquared = maxAxis*maxAxis;
    //如果轴的最大值为-1，则取图像宽高中小的那个一半
    if (maxAxis == -1) {
        if (numRows < numCols) {
            maxAxisSquared = cvRound(0.5*numRows);
        } else {
            maxAxisSquared = cvRound(0.5*numCols);
        }
        maxAxisSquared *= maxAxisSquared;
    }
    
    //定义一些变量 
    int p1, p2, p, x1, y1, x2, y2, x, y;
    double x0, y0, a, b, d, k, dx, dy, fSquared, cosTauSquared, bSquared, alpha;
    //遍历所有非零点 
    for (p1 = 0; p1 < numPixels; p1++) {
    	//取另一个非零点（因为点的顺序不影响，所以，已经处理过的非零点对不再处理） 
        for (p2 = 0; p2 < p1; p2++) {
        	//把点(x1,y1),(x2,y2)取出 
            x1 = nonzeroXIndices[p1];
            y1 = nonzeroYIndices[p1];
            x2 = nonzeroXIndices[p2];
            y2 = nonzeroYIndices[p2];
            
            //计算长半轴 
            dx = x1-x2;
            dy = y1-y2;
            a = 0.5*std::sqrt(dx*dx+dy*dy);
            //看第一个轴的最大值是否为-1，用于判断计算得到的a是否符合要求 
            bool bMaxFirst = (maxFirst == -1);
            if (maxFirst != -1) {
                bMaxFirst = (a <= maxFirst);
            }
            //如果计算到的长半轴a符合要求（大于最小轴的值，并小于第一个轴最大值） 
            if (a > minAxis && bMaxFirst) {
            	//计算中心点 
                x0 = 0.5*(x1+x2);
                y0 = 0.5*(y1+y2);
                //对非零点对，遍历所有非零的点计算短半轴 
                for (p = 0; p < numPixels; p++) {
                	//取出假设的位于椭圆上的点 
                    x = nonzeroXIndices[p];
                    y = nonzeroYIndices[p];
                    dx = x-x0;
                    dy = y-y0;
                    //计算此点到中心点的距离d 
                    d = std::sqrt(dx*dx+dy*dy);
                    //如果此距离大于最小轴的值，则认为此点有可能为椭圆上的点 
                    if (d > minAxis) {
                        dx = x-x2;
                        dy = y-y2;
                        //计算f的值的平方 
                        fSquared = dx*dx+dy*dy;
                        //计算cosTau 
                        cosTauSquared = (a*a+d*d-fSquared) / (2*a*d);
                        //再计算cosTau的平方 
                        cosTauSquared *= cosTauSquared;
                        //计算b平方的分母 
                        k = a*a-d*d*cosTauSquared;
                        //让b要在合理的条件内 
                        if (k > 0 && cosTauSquared < 1) {
                        	//计算b的平方 
                            bSquared = a*a*d*d*(1-cosTauSquared) / k;
                            //计算到的b的平方要小于最大轴的值的平方 
                            if (bSquared <= maxAxisSquared) {
                            	//把符合要求的短半轴的平方值放入累加器列表中 
                                accumulator.push_back(bSquared);
                            }
                        }
                    }
                }
                
                //对非零点对计算完所有短半轴b，看符合要求的值有没有，即累加器列表里是否有值 
                if (accumulator.size() > 0) {
                	//有值，则说明有可能有相应的短半轴与这一对非零点组成的长轴构成一个椭圆 
                    std::vector<double> bins;
                    //下面是对累加器列表中的短半轴平方的数据进行处理，看能不能找到合适的短半轴，
					//首先是利用精度来创建格子 
                    bins = _create_bins(accumulator, binSize);
                    std::vector<int> hist;
                    //然后按创建的格子计算累计器列表的直方图 
                    hist = _histogram(accumulator, bins);
                    //把直方图计数数据中的最大值找出来 
                    std::vector<int>::iterator histMaxIter = std::max_element(hist.begin(), hist.end());
                    int histMax = *histMaxIter;
                    //如果这个最大计数达到了要求的阈值，则说明这个短半轴的计算是可以构成一个椭圆的 
                    if (histMax > threshold) {
                    	//计算倾斜角度 
                        alpha = std::atan2(y2-y1, x2-x1);
                        //把直方图最大值对应的格子索引计算出来 
                        int maxIndex = histMaxIter - hist.begin();
                        //短半轴的值即为最大直方图计数值对应的格子的值的开方 
                        b = std::sqrt(bins[maxIndex]);
                        if (alpha != 0) {
                        	//在OpenCV中，y轴是向下的，把角度改一下（其实问题不大) 
                            alpha = CV_PI - alpha;
                        }
                        //最后保证计算得到的长半轴和短半轴都大于0，则认为确定了一个椭圆 
                        if (a > 0 && b > 0) {
                        	//把确定的椭圆保存 
                            Ellipse ellipse;
                            ellipse.accumulator = histMax;
                            ellipse.xc = x0;
                            ellipse.yc = y0;
                            ellipse.a = a;
                            ellipse.b = b;
                            ellipse.angle = alpha * 180 / CV_PI;
                            results.push_back(ellipse);
                        }
                    }
                    //累加器清零再进行下一对非零点 
                    accumulator.clear();
                }
            }
        }
    }
    
    return results;
}

霍夫椭圆（二）

//直接上函数，传进来源图像的高、宽、轮廓线的点，阈值，然后输出就是椭圆的参数 
static bool _hough_ellipse2(int numRows, int numCols, const std::vector<cv::Point>& contour, 
double threshold, cv::Point& center, float& a, float& b, float& theta) {
    //定义最大距离的图，即源图上每个点到轮廓线（这个轮廓线是先被假定是一个椭圆）的最大距离的图 
	cv::Mat maxDists(numRows, numCols, CV_32FC1, cv::Scalar::all(0));
	//任一点到轮廓线所有点的距离列表 
    std::vector<float> dists(contour.size());
    //图像上每一个点到轮廓线上的点的距离计算 
    for (int i = 0; i < numRows; i++) {
        for (int j = 0; j < numCols; j++) {
            for (int k = 0; k < contour.size(); k++) {
            	//计算图像上的点到轮廓线上的点的距离，并放入列表中 
                float dx = j-contour[k].x;
                float dy = i-contour[k].y;
                float dist = std::sqrt(dx*dx+dy*dy);
                dists[k] = dist;
            }
            //从距离列表中找出最大的那个距离 
            std::vector<float>::iterator iter = std::max_element(dists.begin(), dists.end());
            float maxDist = *iter;
            //把这个最大距离放到图中相应位置 
            maxDists.ptr<float>(i)[j] = maxDist;
        }
    }
    
    //有了最大距离的图，找到这个最大距离图中最小的那个距离，及这个最小距离所在的位置 
    double minVal;
    cv::Point minLoc;
    cv::minMaxLoc(maxDists, &minVal, NULL, &minLoc);
    //根据定理，最大距离中的最小那个对应椭圆的长半轴，而对应的位置就是椭圆中心 
    center = minLoc;
    a = (float)minVal;
    
    //接下来由短半轴和角度定义一个霍夫空间 
    cv::Mat houghSpace(std::floor(a+1), 180, CV_32SC1, cv::Scalar::all(0));
    //因为轮廓线被假定为椭圆了，所以把轮廓线上的点，再加上离散的角度，带入椭圆公式计算短半轴 
    for (int i = 0; i < contour.size(); i++) {
        for (int j = 0; j < 180; j++) {
        	//角度转为弧度 
            double angle = j*CV_PI/180;
            //中心点，前面已经确认 
            double p = center.x;
            double q = center.y;
            //轮廓上的点 
            double x = contour[i].x;
            double y = contour[i].y;
            //计算sin，cos 
            double cosTheta = std::cos(angle);
            double sinTheta = std::sin(angle);
            //公式中加号左边那一项 
            double part1 = std::pow((x-p)*cosTheta+(y-q)*sinTheta, 2) / (a*a);
            //加号右边那一项的分子 
            double part2 = std::pow(-(x-p)*sinTheta+(y-q)*cosTheta, 2);
            //计算短半轴，要取到整数，所以用了floor这个函数，先进行了加1 
            int B = std::floor(std::sqrt(part2 / (1-part1)) + 1);
            //计算得到的短半轴大于0，且不大于长半轴，则霍夫空间对应短半轴和角度位置计算加1 
            if (B > 0 && B <= a) {
                houghSpace.ptr<int>(B)[j] += 1;
            }
        }
    }
    
    //得到霍夫空间的累加计数，下面把累加计数的最大值找出来 
    double maxVal;
    cv::Point maxLoc;
    cv::minMaxLoc(houghSpace, NULL, &maxVal, NULL, &maxLoc);
    //如果最大计数大于给定的阈值，则说明这条轮廓线有成为椭圆的潜质 
    if (maxVal >= threshold) {
    	//短半轴就是最大计数对应位置的y，角度就是x 
        b = maxLoc.y;
        theta = maxLoc.x;
        //最后确保长短半轴都是大于1的，这才说明成功确认一条轮廓线为一个椭圆 
        return (a > 1 && b > 1);
    }
    return false;
}

代码

//网上代码：后期发现不需要也就不改了……
#include <iostream>
#include <algorithm>
#include <opencv2/core.hpp>
#include <opencv2/imgcodecs.hpp>
#include <opencv2/imgproc.hpp>
#include <opencv2/highgui.hpp>

//定义椭圆结构体，成员包含累积计数，中心点，长短半轴，倾斜角度 
typedef struct Ellipse {
    int accumulator;
    double xc, yc;
    double a, b;
    double angle;
} Ellipse;

//创建格子，传累加器列表和精度 
static std::vector<double> _create_bins(const std::vector<double>& acc, double binSize) {
	//把累加器列表中最大值取出，即计算得到的短半轴平方的所有值中最大的值 
    std::vector<double>::const_iterator maxPos = std::max_element(acc.begin(), acc.end());
    double maxVal = *maxPos;
    //要返回的格子列表 
    std::vector<double> v;
    //格子是从0开始，到最大值，间隔为精度 
    for (double d = 0; d < maxVal + binSize; d += binSize) {
        v.push_back(d);
    }
    return v;
}

//然后是计算直方图数据的函数，传入累加器列表和格子列表，返回就是累加器列表中的值落入相应的格子的计算的列表 
static std::vector<int> _histogram(const std::vector<double>& acc, const std::vector<double>& bins) {
	//直方图的数据列表，大小比格子的边界数少1 
    std::vector<int> hist(bins.size()-1);
    //直方图数据初始化为0 
    for (int i = 0; i < hist.size(); i++) {
        hist[i] = 0;
    }
    //遍历累加器列表中的数据，看落入了哪个格子 
    for (int i = 0; i < acc.size(); i++) {
        double a = acc[i];
        //看落入哪个格子 
        for (int j = 0; j < bins.size()-1; j++) {
        	//格子的左边值和右边值，除了最后一个格子，都是左闭右开 
            double left = bins[j];
            double right = bins[j+1];
            //如果值落入相应的格子，相应的直方图计算加1 
            if (a >= left && a < right) {
                hist[j]++;
                break;
            }
        }
        //如果为最后一个格子的边界，最后一个格子的计算也要加1，最后一个格子是左闭右也闭 
        if (a == bins[bins.size()-1]) {
            hist[bins.size()-2]++;
        }
    }
    return hist;
}

//定义算法的实现函数，返回检测到的椭圆列表
//参数：传入的边缘图，
//累计计数的阈值（即达到什么值才认为有效），
//精度（主要用来定义直方图的格子），
//轴的最小值，
//轴的最大值（如果为-1，则取图像宽高中小的那个一半），
//第一个轴的最大值（如果为-1，则无限制） 
static std::vector<Ellipse> _hough_ellipse(const cv::Mat& edges, int threshold, double accuracy,
 int minAxis, int maxAxis=-1, int maxFirst=-1) {
    //定义返回的列表 
	std::vector<Ellipse> results;
    
    //拿到图像的宽高 
    int numRows = edges.rows;
    int numCols = edges.cols;
    //边缘图像的非零点位置 
    std::vector<int> nonzeroXIndices, nonzeroYIndices;
    //把边缘图像的非零点找出来 
    for (int i = 0; i < numRows; i++) {
        for (int j = 0; j < numCols; j++) {
            if (edges.ptr<uchar>(i)[j] > 0) {
                nonzeroXIndices.push_back(j);
                nonzeroYIndices.push_back(i);
            }
        }
    }
    //总共有多少非零点 
    int numPixels = nonzeroXIndices.size();
    //累加器列表 
    std::vector<double> accumulator;
    //精度，即直方图的格子宽度，这里取平方，因为下面的算法计算短半轴时为了减少计算量也是用平方（不开方） 
    double binSize = accuracy*accuracy;
    //最大的半轴值，也取平方以避免开方 
    double maxAxisSquared = maxAxis*maxAxis;
    //如果轴的最大值为-1，则取图像宽高中小的那个一半
    if (maxAxis == -1) {
        if (numRows < numCols) {
            maxAxisSquared = cvRound(0.5*numRows);
        } else {
            maxAxisSquared = cvRound(0.5*numCols);
        }
        maxAxisSquared *= maxAxisSquared;
    }
    
    //定义一些变量 
    int p1, p2, p, x1, y1, x2, y2, x, y;
    double x0, y0, a, b, d, k, dx, dy, fSquared, cosTauSquared, bSquared, alpha;
    //遍历所有非零点 
    for (p1 = 0; p1 < numPixels; p1++) {
    	//取另一个非零点（因为点的顺序不影响，所以，已经处理过的非零点对不再处理） 
        for (p2 = 0; p2 < p1; p2++) {
        	//把点(x1,y1),(x2,y2)取出 
            x1 = nonzeroXIndices[p1];
            y1 = nonzeroYIndices[p1];
            x2 = nonzeroXIndices[p2];
            y2 = nonzeroYIndices[p2];
            
            //计算长半轴 
            dx = x1-x2;
            dy = y1-y2;
            a = 0.5*std::sqrt(dx*dx+dy*dy);
            //看第一个轴的最大值是否为-1，用于判断计算得到的a是否符合要求 
            bool bMaxFirst = (maxFirst == -1);
            if (maxFirst != -1) {
                bMaxFirst = (a <= maxFirst);
            }
            //如果计算到的长半轴a符合要求（大于最小轴的值，并小于第一个轴最大值） 
            if (a > minAxis && bMaxFirst) {
            	//计算中心点 
                x0 = 0.5*(x1+x2);
                y0 = 0.5*(y1+y2);
                //对非零点对，遍历所有非零的点计算短半轴 
                for (p = 0; p < numPixels; p++) {
                	//取出假设的位于椭圆上的点 
                    x = nonzeroXIndices[p];
                    y = nonzeroYIndices[p];
                    dx = x-x0;
                    dy = y-y0;
                    //计算此点到中心点的距离d 
                    d = std::sqrt(dx*dx+dy*dy);
                    //如果此距离大于最小轴的值，则认为此点有可能为椭圆上的点 
                    if (d > minAxis) {
                        dx = x-x2;
                        dy = y-y2;
                        //计算f的值的平方 
                        fSquared = dx*dx+dy*dy;
                        //计算cosTau 
                        cosTauSquared = (a*a+d*d-fSquared) / (2*a*d);
                        //再计算cosTau的平方 
                        cosTauSquared *= cosTauSquared;
                        //计算b平方的分母 
                        k = a*a-d*d*cosTauSquared;
                        //让b要在合理的条件内 
                        if (k > 0 && cosTauSquared < 1) {
                        	//计算b的平方 
                            bSquared = a*a*d*d*(1-cosTauSquared) / k;
                            //计算到的b的平方要小于最大轴的值的平方 
                            if (bSquared <= maxAxisSquared) {
                            	//把符合要求的短半轴的平方值放入累加器列表中 
                                accumulator.push_back(bSquared);
                            }
                        }
                    }
                }
                
                //对非零点对计算完所有短半轴b，看符合要求的值有没有，即累加器列表里是否有值 
                if (accumulator.size() > 0) {
                	//有值，则说明有可能有相应的短半轴与这一对非零点组成的长轴构成一个椭圆 
                    std::vector<double> bins;
                    //下面是对累加器列表中的短半轴平方的数据进行处理，看能不能找到合适的短半轴，
					//首先是利用精度来创建格子 
                    bins = _create_bins(accumulator, binSize);
                    std::vector<int> hist;
                    //然后按创建的格子计算累计器列表的直方图 
                    hist = _histogram(accumulator, bins);
                    //把直方图计数数据中的最大值找出来 
                    std::vector<int>::iterator histMaxIter = std::max_element(hist.begin(), hist.end());
                    int histMax = *histMaxIter;
                    //如果这个最大计数达到了要求的阈值，则说明这个短半轴的计算是可以构成一个椭圆的 
                    if (histMax > threshold) {
                    	//计算倾斜角度 
                        alpha = std::atan2(y2-y1, x2-x1);
                        //把直方图最大值对应的格子索引计算出来 
                        int maxIndex = histMaxIter - hist.begin();
                        //短半轴的值即为最大直方图计数值对应的格子的值的开方 
                        b = std::sqrt(bins[maxIndex]);
                        if (alpha != 0) {
                        	//在OpenCV中，y轴是向下的，把角度改一下（其实问题不大) 
                            alpha = CV_PI - alpha;
                        }
                        //最后保证计算得到的长半轴和短半轴都大于0，则认为确定了一个椭圆 
                        if (a > 0 && b > 0) {
                        	//把确定的椭圆保存 
                            Ellipse ellipse;
                            ellipse.accumulator = histMax;
                            ellipse.xc = x0;
                            ellipse.yc = y0;
                            ellipse.a = a;
                            ellipse.b = b;
                            ellipse.angle = alpha * 180 / CV_PI;
                            results.push_back(ellipse);
                        }
                    }
                    //累加器清零再进行下一对非零点 
                    accumulator.clear();
                }
            }
        }
    }
    
    return results;
}

//直接上函数，传进来源图像的高、宽、轮廓线的点，阈值，然后输出就是椭圆的参数 
static bool _hough_ellipse2(int numRows, int numCols, const std::vector<cv::Point>& contour, 
double threshold, cv::Point& center, float& a, float& b, float& theta) {
    //定义最大距离的图，即源图上每个点到轮廓线（这个轮廓线是先被假定是一个椭圆）的最大距离的图 
	cv::Mat maxDists(numRows, numCols, CV_32FC1, cv::Scalar::all(0));
	//任一点到轮廓线所有点的距离列表 
    std::vector<float> dists(contour.size());
    //图像上每一个点到轮廓线上的点的距离计算 
    for (int i = 0; i < numRows; i++) {
        for (int j = 0; j < numCols; j++) {
            for (int k = 0; k < contour.size(); k++) {
            	//计算图像上的点到轮廓线上的点的距离，并放入列表中 
                float dx = j-contour[k].x;
                float dy = i-contour[k].y;
                float dist = std::sqrt(dx*dx+dy*dy);
                dists[k] = dist;
            }
            //从距离列表中找出最大的那个距离 
            std::vector<float>::iterator iter = std::max_element(dists.begin(), dists.end());
            float maxDist = *iter;
            //把这个最大距离放到图中相应位置 
            maxDists.ptr<float>(i)[j] = maxDist;
        }
    }
    
    //有了最大距离的图，找到这个最大距离图中最小的那个距离，及这个最小距离所在的位置 
    double minVal;
    cv::Point minLoc;
    cv::minMaxLoc(maxDists, &minVal, NULL, &minLoc);
    //根据定理，最大距离中的最小那个对应椭圆的长半轴，而对应的位置就是椭圆中心 
    center = minLoc;
    a = (float)minVal;
    
    //接下来由短半轴和角度定义一个霍夫空间 
    cv::Mat houghSpace(std::floor(a+1), 180, CV_32SC1, cv::Scalar::all(0));
    //因为轮廓线被假定为椭圆了，所以把轮廓线上的点，再加上离散的角度，带入椭圆公式计算短半轴 
    for (int i = 0; i < contour.size(); i++) {
        for (int j = 0; j < 180; j++) {
        	//角度转为弧度 
            double angle = j*CV_PI/180;
            //中心点，前面已经确认 
            double p = center.x;
            double q = center.y;
            //轮廓上的点 
            double x = contour[i].x;
            double y = contour[i].y;
            //计算sin，cos 
            double cosTheta = std::cos(angle);
            double sinTheta = std::sin(angle);
            //公式中加号左边那一项 
            double part1 = std::pow((x-p)*cosTheta+(y-q)*sinTheta, 2) / (a*a);
            //加号右边那一项的分子 
            double part2 = std::pow(-(x-p)*sinTheta+(y-q)*cosTheta, 2);
            //计算短半轴，要取到整数，所以用了floor这个函数，先进行了加1 
            int B = std::floor(std::sqrt(part2 / (1-part1)) + 1);
            //计算得到的短半轴大于0，且不大于长半轴，则霍夫空间对应短半轴和角度位置计算加1 
            if (B > 0 && B <= a) {
                houghSpace.ptr<int>(B)[j] += 1;
            }
        }
    }
    
    //得到霍夫空间的累加计数，下面把累加计数的最大值找出来 
    double maxVal;
    cv::Point maxLoc;
    cv::minMaxLoc(houghSpace, NULL, &maxVal, NULL, &maxLoc);
    //如果最大计数大于给定的阈值，则说明这条轮廓线有成为椭圆的潜质 
    if (maxVal >= threshold) {
    	//短半轴就是最大计数对应位置的y，角度就是x 
        b = maxLoc.y;
        theta = maxLoc.x;
        //最后确保长短半轴都是大于1的，这才说明成功确认一条轮廓线为一个椭圆 
        return (a > 1 && b > 1);
    }
    return false;
}

int main(int argc, char **argv) {
    /*
	//测试代码，读入一张图片 
	cv::Mat src = cv::imread("/home/xlll/Downloads/opencv/samples/data/pic3.png", cv::IMREAD_COLOR);
    if (src.empty()) {
        std::cout << "failed to read image!" << std::endl;
        return EXIT_FAILURE;
    }
    //转为灰度图 
    cv::Mat gray;
    cv::cvtColor(src, gray, cv::COLOR_BGR2GRAY);
    //高斯平滑去一下噪声 
    cv::Mat gb;
    cv::GaussianBlur(gray, gb, cv::Size(7, 7), 2, 2);
    int kernelSize = 3;
    //使用Canny算法检测这张图的边缘图 
    double cannyThresh = 100;
    cv::Mat edges, dx, dy;
    cv::Sobel(gb, dx, CV_16SC1, 1, 0, kernelSize, 1, 0, cv::BORDER_REPLICATE);
    cv::Sobel(gb, dy, CV_16SC1, 0, 1, kernelSize, 1, 0, cv::BORDER_REPLICATE);
    cv::Canny(dx, dy, edges, cannyThresh/2, cannyThresh, false);
    
    //用于霍夫椭圆检测的参数，这个算法的参数不好调，这几个参数是根据这张图调好的值 
    int threshold = 135;
    double accuracy = 15;
    int minAxis = 30;
    int maxFirst = 48;
    int maxAxis = 58;
    //调用算法，返回检测到的椭圆 
    std::vector<Ellipse> ellipses;
    ellipses = _hough_ellipse(edges, threshold, accuracy, minAxis, maxAxis, maxFirst);
    //看检测到几个 
    std::cout << "Found " << ellipses.size() << " ellipses." << std::endl;
    
    //原图克隆一份，用于画椭圆 
    cv::Mat dst = src.clone();
    //把所有检测到的椭圆画到图上，并把椭圆相关参数打印出来 
    for (int i = 0; i < ellipses.size(); i++) {
        Ellipse ellipse = ellipses[i];
        cv::ellipse(dst, cv::Point(cvRound(ellipse.xc), cvRound(ellipse.yc)), cv::Size(cvRound(ellipse.a), cvRound(ellipse.b)), ellipse.angle, 0, 360, cv::Scalar(255, 0, 255), 2, cv::LINE_AA);
        std::cout << "Ellipse " << i << ": xc = " << cvRound(ellipse.xc) << ", yc = " << cvRound(ellipse.yc) << ", a = " << cvRound(ellipse.a) << ", b = " << cvRound(ellipse.b) << ", angle = " << ellipse.angle << ", accumulator = " << ellipse.accumulator << std::endl;
    }
    //显示图像 
    cv::imshow("src", src);
    cv::imshow("edges", edges);
    cv::imshow("dst", dst);
    cv::waitKey(0);*/
    
    
    //测试代码，读入一张图片 
    cv::Mat src = cv::imread("/home/xlll/Downloads/opencv/samples/data/detect_blob.png", cv::IMREAD_COLOR);
    if (src.empty()) {
        std::cout << "failed to read image!" << std::endl;
        return EXIT_FAILURE;
    }
    //转为灰度图 
    cv::Mat gray;
    cv::cvtColor(src, gray, cv::COLOR_BGR2GRAY);   
    cv::Mat gb;
    //高斯平滑一下 
    cv::GaussianBlur(gray, gb, cv::Size(3, 3), 2, 2);
    //使用Canny算法检测这张图的边缘图 
    int kernelSize = 3;
    double cannyThresh = 35;
    cv::Mat edges, dx, dy;
    cv::Sobel(gb, dx, CV_16SC1, 1, 0, kernelSize, 1, 0, cv::BORDER_REPLICATE);
    cv::Sobel(gb, dy, CV_16SC1, 0, 1, kernelSize, 1, 0, cv::BORDER_REPLICATE);
    cv::Canny(dx, dy, edges, cannyThresh/2, cannyThresh, false);
    
    //查找边缘图的轮廓线 
    std::vector<std::vector<cv::Point>> contours;
    cv::findContours(edges, contours, cv::RETR_LIST, cv::CHAIN_APPROX_NONE);
    //克隆图片用于输出 
    cv::Mat dst = src.clone();
    for (int i = 0; i < contours.size(); i++) {
        float a, b, theta;
        cv::Point center;
        //阈值定义，可根据实际情况更改 
        double threshold = contours[i].size() * 0.33;
        //如果成功检测为椭圆，把椭圆画出来，并输出椭圆相应的参数 
        if (_hough_ellipse2(dst.rows, dst.cols, contours[i], threshold, center, a, b, theta)) {
            cv::ellipse(dst, cv::Point(center.x, center.y), cv::Size(cvRound(a), cvRound(b)), theta, 0, 360, cv::Scalar(255, 0, 255), 2, cv::LINE_AA);
            std::cout << "center = (" << center.x << ", " << center.y << "), a = " << a << ", b = " << b << ", theta = " << theta << std::endl;
        }
    }
    //显示相应的图像 
    cv::imshow("src", src);
    cv::imshow("edges", edges);
    cv::imshow("dst", dst);
    cv::waitKey(0);
    
    return 0;
}

本文链接： https://zhang0224gz.github.io/2022/07/04/智能车/霍夫变换检测椭圆/

版权声明： 本博客所有文章除特别声明外，均采用 CC BY 4.0 CN协议许可协议。转载请注明出处！

Eureka

Study Notes & Supplement

霍夫变换检测椭圆

霍夫变换（一）

霍夫椭圆（二）

代码

EurekaStudy Notes & Supplement