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

Question

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

这里有一张图片来说明：

每个表格由 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 的更 analytical/elegant 的解决方案。或者还有什么我可以做的吗？

谢谢！我希望我说清楚了。

Answer 1

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

想象一条从上到下、从左到右移动的线，无论哪个方向对您来说都更容易处理。假设我们将它从上到下移动。

在每个点，与线条重叠的形状都会生成间隔（在线条遇到形状的地方打开，在线条离开形状的地方闭合。给定间隔，您可以轻松计算出读取了多少线条（在任何区间之外）。这是 O(N) 的复杂度。

好的，现在我们需要将它从上到下移动并对面积求和。但是您不想逐个像素地移动它，因为这会使它依赖于 DPI。所以请注意，随着线的移动，间隔会略微移动，但它们不会改变 shape/unless 形状在该点相交，线遇到一个新的形状，而你在上面留下一个形状。

取每个形状的最小和最大 y 值并对它们进行排序。您还需要每对重叠形状之间的交集。但是您不能计算它们，而只能计算扫描线上彼此靠近的那些（将它们放在堆中，因为您需要在每一步中找到下一个）。

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

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

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

Answer 2

这是几何问题。您需要从总面积总和中删除重叠面积。这不是一件容易的事。您可以尝试通过椭圆弧曲线之间的积分来解决它。对于 2 个椭圆就这么简单，但对于更多椭圆，您首先需要决定使用哪个弧以及它是内部还是外部。这可能会变成乏味的代码。

相反，您可以尝试将场景划分为宽度足够小（积分精度）的垂直切片，然后仅使用命中测试 determine/sum 所有未填充区域 光线投射 .

按dx步骤处理场景

dx 步骤将决定结果的准确性。
为每个x收集所有交点

这是简单的 O(N) 搜索，其中 N 是省略号的数量。只需根据当前 x 计算每个椭圆的两个 y 坐标并将其添加到列表（如果有效）。
删除另一个椭圆内的所有交点

这是 O(n^2)，其中 n 是交叉点的数量。只需测试任何交点是否在任何椭圆交点内。如果是，请将其标记为删除……并在完成所有操作后将其删除。不要在搜索过程中删除它们，否则会使结果无效。这样你的列表将只包含外部交点。当心可能存在重复！！！
按y
对路口进行排序
这将加快/简化下一个过程。
从 y=y_start
投射光线
找到路径中的第一个交点

它只是列表中的第一个交点。只需检查它是否在积分范围内，如果不在 skip/stop 内。如果有效，则将距离添加到最后一个点（y_start 或上次迭代的 y）。也不要忘记将 y 值限制为 y_end 否则你可以将屏幕上不可见的区域添加到该区域...
找到这个椭圆的终点

简单地增加指向交点列表的索引，直到 y 坐标跳转（以处理重复值）。如果列表末尾使用 y_end 作为 y 值...添加到区域并停止。
循环 #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 的结果：

WindowCaption中的数字是[pixels^2]中的计算面积，大致是window面的31%。您可以使用任何集成规则或 dx 步骤。您也可以针对不同的事物更改冒泡排序，但通常对于如此小的 n 冒泡会打败其他人......并且代码足够简单。

代码本身尚未优化...

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

Sum area between groups of overlapping forms

python

algorithm

math

image-processing

ellipse