c++ - 如何从椭圆图像中找到上下弧的长度

Question

在这里，我尝试使用图像 vector (图像轮廓)找到上圆弧和下圆弧，但无法给出提取结果。建议使用其他方法从图像及其长度中查找上下弧。

这是我的代码

    Mat image =cv::imread("thinning/20d.jpg");
    int i=0,j=0,k=0,x=320;
    for(int y = 0; y < image.rows; y++)
    {
    if(image.at<Vec3b>(Point(x, y))[0] >= 250 && image.at<Vec3b>(Point(x, y))[1] >= 250 && image.at<Vec3b>(Point(x, y))[2] >= 250){
              qDebug()<<x<<y;
              x1[i]=x;
              y1[i]=y;
              i=i+1;
    }
    }
    for(i=0;i<=1;i++){
      qDebug()<<x1[i]<<y1[i];
    }
    qDebug()<<"UPPER ARC";
    for(int x = 0; x < image.cols; x++)
    {
      for(int y = 0; y <= (y1[0]+20); y++)
      {
          if(image.at<Vec3b>(Point(x, y))[0] >= 240 && image.at<Vec3b>(Point(x, y))[1] >= 240 && image.at<Vec3b>(Point(x, y))[2] >= 240){
              x2[j]=x;
              y2[j]=y;
              j=j+1;
            qDebug()<<x<<y;
          }}
    }   
    qDebug()<<"Lower ARC";
    for(int x = 0; x < image.cols; x++)
    {
      for(int y = (y1[1]-20); y <= image.rows; y++)
      {
          if(image.at<Vec3b>(Point(x, y))[0] >= 240 && image.at<Vec3b>(Point(x, y))[1] >= 240 && image.at<Vec3b>(Point(x, y))[2] >= 240){
              x3[k]=x;
              y3[k]=y;
              k=k+1;
             qDebug()<<x<<y;
          }}
   }

通过上述代码，我得到了坐标，通过使用坐标点，我可以找到圆弧的长度，但是它与提取结果不匹配。

这是实际图像:

图片1:



经过细化后，我得到:



预期产量:

Answer 1

由于您无法定义确切的上/下弧线，因此我假设您通过穿过椭圆的中点的水平线将椭圆切成两半。如果不是这种情况，那么您必须自己调整它...好吧，现在该怎么做:

将图像二值化

当您提供JPG时，颜色会失真，因此只有黑白

将边框缩小到1像素

用白色填充内部，然后将不与任何黑色像素相邻的所有白色像素重新着色为一些未使用或黑色。还有许多其他变种如何实现这一目标...

找到边界框

搜索所有像素，并记住所有白色像素的最小，最大x,y坐标。我们称它们为x0,y0,x1,y1。

椭圆计算中心

只需找到边界框的中间点

cx=(x0+x1)/2
cy=(y0+y1)/2

计算每个椭圆弧的像素

对每个弧都有一个计数器，并为具有y<=cy的任何白色像素简单地增加高弧计数器，如果是y>=cy则降低它。如果坐标系不同，则条件可以相反。

查找椭圆参数

只需找到最接近(cx,cy)的白色像素，这将是次半轴b的端点，就称它为(bx,by)。还要找到(cx,cy)最远的白色像素，它将成为主要的半轴端点(ax,ay)。它们与中心之间的距离将为您提供a,b，而中心减去的位置将为您提供椭圆旋转的 vector 。可以通过atan2获得 Angular ，也可以像我一样使用基 vector 。您可以按点积测试正交性。最接近和最远的点可能超过2点。在这种情况下，您应该找到每个组的中间以提高精度。

集成拟合的椭圆

您首先需要使用y=cy查找椭圆点所在的 Angular ，然后对这两个 Angular 之间的椭圆进行积分。另一半是相同的，只是集成了angles + PI。要确定哪一半只是 Angular 范围中间的计算点，并根据y>=cy ...

进行确定

[Edit2]为此，我为此更新了C++代码:

    picture pic0,pic1,pic2;
        // pic0 - source
        // pic1 - output
    float a,b,a0,a1,da,xx0,xx1,yy0,yy1,ll0,ll1;
    int x,y,i,threshold=127,x0,y0,x1,y1,cx,cy,ax,ay,bx,by,aa,bb,dd,l0,l1;
    pic1=pic0;
    // bbox,center,recolor (white,black)
    x0=pic1.xs; x1=0;
    y0=pic1.ys; y1=0;
    for (y=0;y<pic1.ys;y++)
     for (x=0;x<pic1.xs;x++)
      if (pic1.p[y][x].db[0]>=threshold)
        {
        if (x0>x) x0=x;
        if (y0>y) y0=y;
        if (x1<x) x1=x;
        if (y1<y) y1=y;
        pic1.p[y][x].dd=0x00FFFFFF;
        } else pic1.p[y][x].dd=0x00000000;
    cx=(x0+x1)/2; cy=(y0+y1)/2;
    // fill inside (gray) left single pixel width border (thining)
    for (y=y0;y<=y1;y++)
        {
        for (x=x0;x<=x1;x++) if (pic1.p[y][x].dd)
            {
            for (i=x1;i>=x;i--) if (pic1.p[y][i].dd)
                {
                for (x++;x<i;x++) pic1.p[y][x].dd=0x00202020;
                break;
                }
            break;
            }
        }
    for (x=x0;x<=x1;x++)
        {
        for (y=y0;y<=y1;y++) if (pic1.p[y][x].dd) { pic1.p[y][x].dd=0x00FFFFFF; break; }
        for (y=y1;y>=y0;y--) if (pic1.p[y][x].dd) { pic1.p[y][x].dd=0x00FFFFFF; break; }
        }
    // find min,max radius (periaxes)
    bb=pic1.xs+pic1.ys; bb*=bb; aa=0;
    ax=cx; ay=cy; bx=cx; by=cy;
    for (y=y0;y<=y1;y++)
     for (x=x0;x<=x1;x++)
      if (pic1.p[y][x].dd==0x00FFFFFF)
        {
        dd=((x-cx)*(x-cx))+((y-cy)*(y-cy));
        if (aa<dd) { ax=x; ay=y; aa=dd; }
        if (bb>dd) { bx=x; by=y; bb=dd; }
        }
    aa=sqrt(aa); ax-=cx; ay-=cy;
    bb=sqrt(bb); bx-=cx; by-=cy;
    //a=float((ax*bx)+(ay*by))/float(aa*bb);    // if (fabs(a)>zero_threshold) not perpendicular semiaxes

    // separate/count upper,lower arc by horizontal line
    l0=0; l1=0;
    for (y=y0;y<=y1;y++)
     for (x=x0;x<=x1;x++)
      if (pic1.p[y][x].dd==0x00FFFFFF)
        {
        if (y>=cy) { l0++; pic1.p[y][x].dd=0x000000FF; } // red
        if (y<=cy) { l1++; pic1.p[y][x].dd=0x00FF0000; } // blue
        }
    // here is just VCL/GDI info layer output so you can ignore it...

    // arc separator axis
    pic1.bmp->Canvas->Pen->Color=0x00808080;
    pic1.bmp->Canvas->MoveTo(x0,cy);
    pic1.bmp->Canvas->LineTo(x1,cy);

    // draw analytical ellipse to compare
    pic1.bmp->Canvas->Pen->Color=0x0000FF00;
    pic1.bmp->Canvas->MoveTo(cx,cy);
    pic1.bmp->Canvas->LineTo(cx+ax,cy+ay);
    pic1.bmp->Canvas->MoveTo(cx,cy);
    pic1.bmp->Canvas->LineTo(cx+bx,cy+by);
    pic1.bmp->Canvas->Pen->Color=0x00FFFF00;
    da=0.01*M_PI;   // dash step [rad]
    a0=0.0;         // start
    a1=2.0*M_PI;    // end
    for (i=1,a=a0;i;)
        {
        a+=da; if (a>=a1) { a=a1; i=0; }
        x=cx+(ax*cos(a))+(bx*sin(a));
        y=cy+(ay*cos(a))+(by*sin(a));
        pic1.bmp->Canvas->MoveTo(x,y);
        a+=da; if (a>=a1) { a=a1; i=0; }
        x=cx+(ax*cos(a))+(bx*sin(a));
        y=cy+(ay*cos(a))+(by*sin(a));
        pic1.bmp->Canvas->LineTo(x,y);
        }

    // integrate the arclengths from fitted ellipse
    da=0.001*M_PI;      // integration step [rad] (accuracy)
    // find start-end angles
    ll0=M_PI; ll1=M_PI;
    for (i=1,a=0.0;i;)
        {
        a+=da; if (a>=2.0*M_PI) { a=0.0; i=0; }
        xx1=(ax*cos(a))+(bx*sin(a));
        yy1=(ay*cos(a))+(by*sin(a));
        b=atan2(yy1,xx1);
        xx0=fabs(b-0.0); if (xx0>M_PI) xx0=2.0*M_PI-xx0;
        xx1=fabs(b-M_PI);if (xx1>M_PI) xx1=2.0*M_PI-xx1;
        if (ll0>xx0) { ll0=xx0; a0=a; }
        if (ll1>xx1) { ll1=xx1; a1=a; }
        }
    // [upper half]
    ll0=0.0;
    xx0=cx+(ax*cos(a0))+(bx*sin(a0));
    yy0=cy+(ay*cos(a0))+(by*sin(a0));
    for (i=1,a=a0;i;)
        {
        a+=da; if (a>=a1) { a=a1; i=0; }
        xx1=cx+(ax*cos(a))+(bx*sin(a));
        yy1=cy+(ay*cos(a))+(by*sin(a));
        // sum arc-line sizes
        xx0-=xx1; xx0*=xx0;
        yy0-=yy1; yy0*=yy0;
        ll0+=sqrt(xx0+yy0);
//      pic1.p[int(yy1)][int(xx1)].dd=0x0000FF00; // recolor for visualy check the right arc selection
        xx0=xx1; yy0=yy1;
        }
    // lower half
    a0+=M_PI; a1+=M_PI; ll1=0.0;
    xx0=cx+(ax*cos(a0))+(bx*sin(a0));
    yy0=cy+(ay*cos(a0))+(by*sin(a0));
    for (i=1,a=a0;i;)
        {
        a+=da; if (a>=a1) { a=a1; i=0; }
        xx1=cx+(ax*cos(a))+(bx*sin(a));
        yy1=cy+(ay*cos(a))+(by*sin(a));
        // sum arc-line sizes
        xx0-=xx1; xx0*=xx0;
        yy0-=yy1; yy0*=yy0;
        ll1+=sqrt(xx0+yy0);
//      pic1.p[int(yy1)][int(xx1)].dd=0x00FF00FF; // recolor for visualy check the right arc selection
        xx0=xx1; yy0=yy1;
        }
    // handle if the upper/lower parts are swapped
    a=a0+0.5*(a1-a0);
    if ((ay*cos(a))+(by*sin(a))<0.0) { a=ll0; ll0=ll1; ll1=a; }
    // info texts
    pic1.bmp->Canvas->Font->Color=0x00FFFF00;
    pic1.bmp->Canvas->Brush->Style=bsClear;
    x=5; y=5; i=16; y-=i;
    pic1.bmp->Canvas->TextOutA(x,y+=i,AnsiString().sprintf("center = (%i,%i) px",cx,cy));
    pic1.bmp->Canvas->TextOutA(x,y+=i,AnsiString().sprintf("a = %i px",aa));
    pic1.bmp->Canvas->TextOutA(x,y+=i,AnsiString().sprintf("b = %i px",bb));
    pic1.bmp->Canvas->TextOutA(x,y+=i,AnsiString().sprintf("upper = %i px",l0));
    pic1.bmp->Canvas->TextOutA(x,y+=i,AnsiString().sprintf("lower = %i px",l1));
    pic1.bmp->Canvas->TextOutA(x,y+=i,AnsiString().sprintf("upper`= %.3lf px",ll0));
    pic1.bmp->Canvas->TextOutA(x,y+=i,AnsiString().sprintf("lower`= %.3lf px",ll1));
    pic1.bmp->Canvas->Brush->Style=bsSolid;

它对成员使用我自己的图片类:

xs,ys图像

的分辨率

p[y][x].dd像素访问为32bit无符号整数，颜色为

p[y][x].db[4]像素访问作为4 * 8bit无符号整数作为颜色通道

您可以将picture::p成员视为简单的2D数组

union color
    {
    DWORD dd; WORD dw[2]; byte db[4];
    int i; short int ii[2];
    color(){}; color(color& a){ *this=a; }; ~color(){}; color* operator = (const color *a) { dd=a->dd; return this; }; /*color* operator = (const color &a) { ...copy... return this; };*/
    };
int xs,ys;
color p[ys][xs];
Graphics::TBitmap *bmp; // VCL GDI Bitmap object you do not need this...

现在可以不确定每个单元格是否可以作为32位像素p[][].dd或0xAABBGGRR或0xAARRGGBB来访问。您也可以使用p [] []。db [4]作为8位BYTE直接访问通道。
bmp成员是 GDI 位图，因此bmp->Canvas->可访问所有对您不重要的 GDI 东西。

第二张图片的结果如下:

灰色水平线是穿过中心

的圆弧边界线

红色，蓝色是弧形两半(计数时重新着色)

绿色是半轴基础 vector

水色破折号是分析椭圆叠加层，用于比较拟合。

如您所见，拟合度非常接近(+/- 1像素)。所计算的弧长 upper,lower非常接近近似的平均圆半周长(周长)。

您应该添加 a0范围检查来确定起点是上半部还是下半部，因为不会隔离将在主轴的哪一侧找到。这两个半部分的积分几乎相同，并且围绕积分步 0.001*M_PI和每个弧长的 307.3 pixels饱和，这仅是 17和 22像素与直接像素数的差异，这甚至更好，由于混叠，我预计...

对于更偏心的椭圆，拟合度不佳，但结果仍然足够好: