如何从 OpenCV 的 fitEllipse 函数中获取椭圆系数? [英] How can I get ellipse coefficient from fitEllipse function of OpenCV?

查看:35
本文介绍了如何从 OpenCV 的 fitEllipse 函数中获取椭圆系数?的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我想从一张图片中提取红球并得到图片中检测到的椭圆矩阵.

这是我的例子:

我对图片进行阈值处理,使用 findContour() 函数找到红球的轮廓,并使用 fitEllipse() 拟合椭圆.

但我想要的是得到这个椭圆的系数.因为fitEllipse()返回的是一个旋转矩形(RotatedRect),所以需要重新编写这个函数.

一个椭圆可以表示为Ax^2 + By^2 + Cxy + Dx + Ey + F = 0;所以我想得到 u=(A,B,C,D,E,F) 或 u=(A,B,C,D,E) 如果 F 是 1(构造一个椭圆矩阵).

我看了fitEllipse()的源码,一共有三个SVD过程,我想我可以从这三个SVD过程的结果中得到上面的系数.但我很困惑每个 SVD 过程的每个结果(变量 cv::Mat x)代表什么,为什么这里有三个 SVD?

这是这个函数:

cv::RotatedRect cv::fitEllipse( InputArray _points ){垫点 = _points.getMat();int i, n = points.checkVector(2);int depth = points.depth();CV_Assert(n >= 0 && (depth == CV_32F || depth == CV_32S));旋转矩形框;如果( n <5 )CV_Error( CV_StsBadSize, "应该至少有 5 个点来拟合椭圆");//新的拟合椭圆算法,由 Daniel Weiss 博士贡献Point2f c(0,0);双 gfp[5], rp[5], t;const double min_eps = 1e-8;bool is_float = depth == CV_32F;const Point* ptsi = points.ptr();const Point2f* ptsf = points.ptr();AutoBuffer_广告(n*5), _bd(n);double *Ad = _Ad, *bd = _bd;//首先适合参数 A - E垫 A( n, 5, CV_64F, Ad );垫 b( n, 1, CV_64F, bd );垫 x( 5, 1, CV_64F, gfp );for( i = 0; i  min_eps)t = gfp[2]/sin(-2.0 * rp[4]);else//椭圆旋转了 pi/2 的整数倍t = gfp[1] - gfp[0];rp[2] = fabs(gfp[0] + gfp[1] - t);如果( rp[2] > min_eps )rp[2] = std::sqrt(2.0/rp[2]);rp[3] = fabs(gfp[0] + gfp[1] + t);如果( rp[3] > min_eps )rp[3] = std::sqrt(2.0/rp[3]);box.center.x = (float)rp[0] + c.x;box.center.y = (float)rp[1] + c.y;box.size.width = (float)(rp[2]*2);box.size.height = (float)(rp[3]*2);if( box.size.width > box.size.height ){浮动 tmp;CV_SWAP( box.size.width, box.size.height, tmp );box.angle = (float)(90 + rp[4]*180/CV_PI);}如果(框.角度 <-180 )box.angle += 360;if(box.angle > 360 )box.angle -= 360;退货箱;}

源码链接:

因此,您可以从RotatedRect 中获取您需要的参数,而无需更改函数fitEllipse.solve 函数用于求解线性系统或最小二乘法问题.使用 SVD 分解方法,系统可以被过度定义和/或矩阵 src1 可以是奇异的.

有关算法的更多详细信息,您可以查看Fitzgibbon 的论文 提出了这种拟合椭圆方法.

I want to extract the red ball from one picture and get the detected ellipse matrix in picture.

Here is my example:

I threshold the picture, find the contour of red ball by using findContour() function and use fitEllipse() to fit an ellipse.

But what I want is to get coefficient of this ellipse. Because the fitEllipse() return a rotation rectangle (RotatedRect), so I need to re-write this function.

One Ellipse can be expressed as Ax^2 + By^2 + Cxy + Dx + Ey + F = 0; So I want to get u=(A,B,C,D,E,F) or u=(A,B,C,D,E) if F is 1 (to construct an ellipse matrix).

I read the source code of fitEllipse(), there are totally three SVD process, I think I can get the above coefficients from the results of those three SVD process. But I am quite confused what does each result (variable cv::Mat x) of each SVD process represent and why there are three SVD here?

Here is this function:

cv::RotatedRect cv::fitEllipse( InputArray _points )
{
   Mat points = _points.getMat();
   int i, n = points.checkVector(2);
   int depth = points.depth();
   CV_Assert( n >= 0 && (depth == CV_32F || depth == CV_32S));

   RotatedRect box;

   if( n < 5 )
        CV_Error( CV_StsBadSize, "There should be at least 5 points to fit the ellipse" );

    // New fitellipse algorithm, contributed by Dr. Daniel Weiss
    Point2f c(0,0);
    double gfp[5], rp[5], t;
    const double min_eps = 1e-8;
    bool is_float = depth == CV_32F;
    const Point* ptsi = points.ptr<Point>();
    const Point2f* ptsf = points.ptr<Point2f>();

    AutoBuffer<double> _Ad(n*5), _bd(n);
    double *Ad = _Ad, *bd = _bd;

    // first fit for parameters A - E
    Mat A( n, 5, CV_64F, Ad );
    Mat b( n, 1, CV_64F, bd );
    Mat x( 5, 1, CV_64F, gfp );

    for( i = 0; i < n; i++ )
    {
        Point2f p = is_float ? ptsf[i] : Point2f((float)ptsi[i].x, (float)ptsi[i].y);
        c += p;
    }
    c.x /= n;
    c.y /= n;

    for( i = 0; i < n; i++ )
    {
        Point2f p = is_float ? ptsf[i] : Point2f((float)ptsi[i].x, (float)ptsi[i].y);
        p -= c;

        bd[i] = 10000.0; // 1.0?
        Ad[i*5] = -(double)p.x * p.x; // A - C signs inverted as proposed by APP
        Ad[i*5 + 1] = -(double)p.y * p.y;
        Ad[i*5 + 2] = -(double)p.x * p.y;
        Ad[i*5 + 3] = p.x;
        Ad[i*5 + 4] = p.y;
    }

    solve(A, b, x, DECOMP_SVD);

    // now use general-form parameters A - E to find the ellipse center:
    // differentiate general form wrt x/y to get two equations for cx and cy
    A = Mat( 2, 2, CV_64F, Ad );
    b = Mat( 2, 1, CV_64F, bd );
    x = Mat( 2, 1, CV_64F, rp );
    Ad[0] = 2 * gfp[0];
    Ad[1] = Ad[2] = gfp[2];
    Ad[3] = 2 * gfp[1];
    bd[0] = gfp[3];
    bd[1] = gfp[4];
    solve( A, b, x, DECOMP_SVD );

    // re-fit for parameters A - C with those center coordinates
    A = Mat( n, 3, CV_64F, Ad );
    b = Mat( n, 1, CV_64F, bd );
    x = Mat( 3, 1, CV_64F, gfp );
    for( i = 0; i < n; i++ )
    {
        Point2f p = is_float ? ptsf[i] : Point2f((float)ptsi[i].x, (float)ptsi[i].y);
        p -= c;
        bd[i] = 1.0;
        Ad[i * 3] = (p.x - rp[0]) * (p.x - rp[0]);
        Ad[i * 3 + 1] = (p.y - rp[1]) * (p.y - rp[1]);
        Ad[i * 3 + 2] = (p.x - rp[0]) * (p.y - rp[1]);
    }
    solve(A, b, x, DECOMP_SVD);

    // store angle and radii
    rp[4] = -0.5 * atan2(gfp[2], gfp[1] - gfp[0]); // convert from APP angle usage
    if( fabs(gfp[2]) > min_eps )
        t = gfp[2]/sin(-2.0 * rp[4]);
    else // ellipse is rotated by an integer multiple of pi/2
        t = gfp[1] - gfp[0];
    rp[2] = fabs(gfp[0] + gfp[1] - t);
    if( rp[2] > min_eps )
        rp[2] = std::sqrt(2.0 / rp[2]);
    rp[3] = fabs(gfp[0] + gfp[1] + t);
    if( rp[3] > min_eps )
        rp[3] = std::sqrt(2.0 / rp[3]);

    box.center.x = (float)rp[0] + c.x;
    box.center.y = (float)rp[1] + c.y;
    box.size.width = (float)(rp[2]*2);
    box.size.height = (float)(rp[3]*2);
    if( box.size.width > box.size.height )
    {
        float tmp;
        CV_SWAP( box.size.width, box.size.height, tmp );
        box.angle = (float)(90 + rp[4]*180/CV_PI);
    }
    if( box.angle < -180 )
        box.angle += 360;
    if( box.angle > 360 )
        box.angle -= 360;

    return box;
}

The source code link: https://github.com/Itseez/opencv/blob/master/modules/imgproc/src/shapedescr.cpp

解决方案

The function fitEllipse returns a RotatedRect that contains all the parameters of the ellipse.

An ellipse is defined by 5 parameters:

  • xc : x coordinate of the center
  • yc : y coordinate of the center
  • a : major semi-axis
  • b : minor semi-axis
  • theta : rotation angle

You can obtain these parameters like:

RotatedRect e = fitEllipse(points);

float xc    = e.center.x;
float yc    = e.center.y;
float a     = e.size.width  / 2;    // width >= height
float b     = e.size.height / 2;
float theta = e.angle;              // in degrees

You can draw an ellipse with the function ellipse using the RotatedRect:

ellipse(image, e, Scalar(0,255,0)); 

or, equivalently using the ellipse parameters:

ellipse(res, Point(xc, yc), Size(a, b), theta, 0.0, 360.0, Scalar(0,255,0));

If you need the values of the coefficients of the implicit equation, you can do like (from Wikipedia):

So, you can get the parameters you need from the RotatedRect, and you don't need to change the function fitEllipse. The solve function is used to solve linear systems or least-squares problems. Using the SVD decomposition method the system can be over-defined and/or the matrix src1 can be singular.

For more details on the algorithm, you can see the paper of Fitzgibbon that proposed this fit ellipse method.

这篇关于如何从 OpenCV 的 fitEllipse 函数中获取椭圆系数?的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆