其他分享
首页 > 其他分享> > c-如何从椭圆图像中找到上下弧的长度

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

作者:互联网

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

这是我的代码

    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;
          }}
   }

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

这是实际图像:

此搜索:

enter image description here

经过细化后,我得到:

enter image description here

预期产量:

enter image description here

解决方法:

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

>将图像二值化

当您提供JPG时,颜色会失真,因此只有黑白
>将边框缩小到1像素

用白色填充内部,然后将不与任何黑色像素相邻的所有白色像素重新着色为一些未使用或黑色.还有许多其他方式可以实现这一目标…
>找到边界框

搜索所有像素,并记住所有白色像素的min,max 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,而中心减去的位置将为您提供椭圆旋转的矢量.可以通过atan2获得角度,或者像我一样使用基向量.您可以按点积测试正交性.最接近和最远的点可能超过2点.在这种情况下,您应该找到每个组的中间以提高精度.
>整合拟合椭圆

您首先需要找到椭圆点在y = cy处的角度,然后对这两个角度之间的椭圆进行积分.另一半相同,只是积分角PI.为了确定哪一半只是计算角度范围中间的点,并根据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...

其中每个单元可以作为0xAABBGGRR或0xAARRGGBB的32位像素p [] [].dd访问,现在不确定是哪个.您也可以使用p [] [].db [4]作为8位BYTE直接访问通道.

bmp成员是GDI位图,因此bmp-> Canvas->访问所有对您不重要的GDI内容.

这是第二张图片的结果:

example

>灰色水平线是穿过中心的圆弧边界线
>红色,蓝色是弧形两半(计数时重新着色)
> Green是半轴基础向量
> Aqua破折号是分析椭圆叠加层,用于比较拟合度.

如您所见,拟合度非常接近(/ -1像素).计算的上,下圆弧长度非常接近近似的平均圆半周长(周长).

您应该添加一个0范围检查来确定起点是上半部还是下半部,因为没有隔离区可以找到主轴的哪一侧.两半的积分几乎相同,并且在积分步骤0.001 * M_PI左右达到饱和,每弧长约307.3像素,与直接像素数相比仅相差17和22像素,由于混叠,这甚至比我预期的还要好.

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

example

标签:c,opencv,qt,image-processing,ellipse
来源: https://codeday.me/bug/20191013/1908014.html