重叠形式组之间的总和面积

问题描述 投票:3回答:2

我无法找到一种干净的方法来计算给定区域中Python绘制的曲线组之间的剩余区域。

这是一张图片来说明:

Area in Red

每个形式由2个半椭圆组成,具有不同的参数高度参数方程:

X     =Xc+A *cos(Theta)
Y_down=Yc+B1*sin(Theta)
Y_up  =Yc+B2*sin(Theta)

沿着1行(X方向),除了Xc之外,参数方程是相同的。沿Y轴(垂直方向),A随Xc和Yc变化。

所有形式都是通过X轴上的迭代和Y轴上的迭代来完成的。我在绘图中使用Zorder,使它们按照创建的顺序重叠。

问题是,即使我可以计算每个表格的面积,我也看不出如何找到红色区域,因为这些表格在各方面都是重叠的。目前,我可以通过绘制所有曲线并将输出图形和求和二值化来访问红色区域。但我想找到一个更加分析/优雅的解决方案,它不依赖于输出数据的DPI。或者我还能做些什么吗?

谢谢!我希望我很清楚。

python algorithm math image-processing ellipse
2个回答
0
投票

这是几何问题。您需要从总面积总和中删除重叠区域。这不是一件容易的事。您可以尝试通过椭圆弧曲线之间的积分来解决它。对于2个椭圆是很容易的,但是你需要首先确定使用哪个弧,以及它是内部还是外部。这可能会导致繁琐的代码。

相反,您可以尝试将场景划分为具有足够小宽度(整合精度)的垂直切片,并且仅使用命中测试通过光线投射确定/求和所有未填充区域。

img

  1. 一些dx步骤的过程场景 dx步骤将确定结果的准确性。
  2. 为每个x收集所有交叉点 这是简单的O(N)搜索,其中N是椭圆的数量。只需根据当前的y计算每个椭圆的两个x坐标,并将其添加到列表中(如果有效)。
  3. 删除另一个椭圆内的所有交叉点 这是O(n^2),其中n是交叉点的数量。只需测试任何交叉点是否在任何椭圆交叉点内。如果是,请将其标记为删除...并在完成所有操作后删除它们。在搜索过程中不要删除它们,因为这会使结果无效。这样你的列表将只包含外部交叉点。当心可能存在重复!!!
  4. y排序交叉点 这将加速/简化下一个过程。
  5. 来自y=y_start的射线
  6. 在路径中找到第一个交叉点 它只是列表中的第一个交叉点。只需检查它是否在集成范围内,如果没有则跳过/停止。如果有效,则将距离添加到最后一个点(y_start或y,从上次迭代开始)。另外不要忘记将y值钳制到你的y_end,否则你可以在屏幕上添加区域不可见的区域......
  7. 找到这个椭圆的末端 只需增加指向交叉点列表的索引,直到y坐标跳转(处理重复值)。如果列表末尾使用y_end作为y值...添加到区域并停止。
  8. 循环#6直到y_end坐标被击中或交叉

这是我的C ++ / VCL实现:

//---------------------------------------------------------------------------
struct _ellipse
    {
    double x,y,rx,ry0,ry1;
    _ellipse()  {}
    _ellipse(_ellipse& a)   { *this=a; }
    ~_ellipse() {}
    _ellipse* operator = (const _ellipse *a) { *this=*a; return this; }
    //_ellipse* operator = (const _ellipse &a) { ...copy... return this; }
    };
struct _hit
    {
    double y;   // hit y coordinate
    int ix;     // ellipse index
    int f;      // 0 means ry0, 1 means ry1 edge, >=2 means deleted
    _hit()  {}
    _hit(_hit& a) { *this=a; }
    ~_hit() {}
    _hit* operator = (const _hit *a) { *this=*a; return this; }
    //_hit* operator = (const _hit &a) { ...copy... return this; }
    };
const int N=50;
_ellipse ell[N];
//---------------------------------------------------------------------------
void sort_asc_bubble(_hit *a,int n)
    {
    int i,e; _hit q;
    for (e=1;e;n--)
     for (e=0,i=1;i<n;i++)
      if (a[i-1].y>a[i].y)
       { q=a[i-1]; a[i-1]=a[i]; a[i]=q; e=1; }
    }
//---------------------------------------------------------------------------
void init()
    {
    int i; double xs=512,ys=512;
    _ellipse a;
    RandSeed=587654321;
    for (i=0;i<N;i++)
        {
        a.x=xs*Random();
        a.y=ys*Random();
        a.rx =xs*(0.02+(0.25*Random()));
        a.ry0=ys*(0.02+(0.09*Random()));
        a.ry1=ys*(0.02+(0.09*Random()));
        ell[i]=a;
        }
    }
//---------------------------------------------------------------------------
double area()
    {
    double area,dx=1.0;
    double x,y,xx,y0,y1,rxx,ryy0,ryy1;
    _ellipse *p;
    _hit  hit[N+N];     // intersection points
    int  i,j,n;         // n = number of intersections
    int _in;
    for (area=0.0,x=0.0;x<scr.xs;x+=dx)     // all vertical slices
        {
        // [prepare data]
        // O(N) precompute intersection points for ray/ellipses
        for (n=0,p=ell,i=0;i<N;i++,p++)
         if ((x>=p->x-p->rx)&&(x<=p->x+p->rx))
            {
            xx=x-p->x; xx*=xx;
            rxx =p->rx *p->rx ;
            ryy0=p->ry0*p->ry0;
            ryy1=p->ry1*p->ry1;
            hit[n].ix=i; hit[n].f=0; hit[n].y=p->y-sqrt(ryy0*(1.0-(xx/rxx))); n++;
            hit[n].ix=i; hit[n].f=1; hit[n].y=p->y+sqrt(ryy1*(1.0-(xx/rxx))); n++;
            }
        // O(n^2) flag inside edges
        for (i=0;i+1<n;i+=2)
            {
            y0=hit[i+0].y;
            y1=hit[i+1].y;
            for (j=0;j<n;j++)
             if ((i!=j)&&(i+1!=j))
              if ((hit[j].y>y0)&&(hit[j].y<y1))
               hit[j].f|=2;
            }
        // O(n) delete flagged edges
        for (i=0,j=0;i<n;i++)
         if (hit[i].f<2)
          { hit[j]=hit[i]; j++; } n=j;
        // O(n^2) sort y asc and original indexes are in i0,i1
        sort_asc_bubble(hit,n);

        // [integrate area]
        for (i=-1,y0=0.0,y1=0.0;;)
            {
            i++; if (i<n) y=hit[i].y; else y=scr.ys;
            if (y>scr.ys) y=y>scr.ys;
            if (y>y1)
                {
                y0=y1; y1=y;
                area+=y1-y0;
                // debug draw for visual check (render red pixel at x,y0 ... x,y1)
                //for (y=y0;y<=y1;y++) scr.pyx[int(y)][int(x)]=0x00FF0000;
                }
            if (i>=n) break;
            y=hit[i].y;
            for (;i<n;i++) if (hit[i].y>y) { y1=hit[i].y; break; }
            }
        } area*=dx; // rectangular rule
    return area;
    }
//---------------------------------------------------------------------------

这里结果为dx=1.0和分辨率512x512

result

Window Caption中的数字是[pixels^2]中窗口表面的大致31%的计算区域。您可以使用任何集成规则或dx步骤。你也可以改变不同的泡泡排序,但通常对于这样小的n泡泡击败其他人......并且很容易编码。

代码本身尚未优化......


1
投票

这听起来像是线扫描问题。

想象一条线从上到下,从左到右移动,无论方向更容易处理。让我们说我们从上到下移动它。

在每个点上,与线重叠的形状将产生间隔(线在线遇到形状时打开,并在线离开形状的地方闭合。给定间隔,您可以轻松计算读取线的数量(在任何间隔之外)这是O(N)的复杂性。

好的,现在我们需要从上到下移动并总结区域。但是你不想逐个像素地移动它,因为这会使它依赖于DPI。因此请注意,当线移动时,间隔会略微移动,但它们不会改变形状/除非形状在该点相交,否则线会遇到新的形状,并在上面留下形状。

获取每个形状的最小值和最大值并对其进行排序。您还需要在每对重叠形状之间进行交叉。但是,除了计算它们之外,您只能计算扫描线上彼此靠近的那些(将它们保存在堆中,因为您需要在每个步骤中找到下一个)。

因此,一种方法是将线移动一点(到下一个感兴趣的点),计算你刚用线扫过的区域并将其添加到总数中。这将是每行移动O(N ^ 3),O(N)并且您可能有O(N ^ 2)行移动(假设每个形状与每个其他形状重叠)。

如果加快面积计算,你可以更快一点O(N ^ 2 log N)。您需要计算的唯一部分是交换的间隔周围。这意味着,您需要记住的每个空闲间隔是您上次更新的。

我把数学留给你,但这很简单。如果你考虑一个区间,它基本上是一个由两个椭圆填充的矩形。矩形的面积很简单,你只需要在椭圆切片之外添加它们。

© www.soinside.com 2019 - 2024. All rights reserved.