如何找到覆盖多边形的圆的最小扇区?

问题描述

给定圆C:(O,r)和多边形P,我如何找到覆盖P的C的最小扇区?

enter image description here

假设圆的半径足够大,所以主要问题是找到扇形的初始和最终角度。

我试图从圆心向多边形的每个角度绘制射线,并检查射线是否与多边形重叠。但是可能有两条以上的光线仅接触多边形。由于双精度,我不能依靠基于它们的方向角选择唯一的射线。因此,在接触光线列表中找到最小和最大角度是没有用的。除此之外,我在选择由两个终端角度创建的扇区中遇到任何问题,因为当由atan2计算时,初始角度可能大于最终角度。

enter image description here

那么找到这样一个部门的正确方法是什么?

编辑: 三个示例多边形(WKT格式):

polyGON((52.87404 30.85613,42.55699 28.46292,41.54373 24.319989,53.57623 21.300564,62.94891 28.46292,49.39652 27.550071,52.87404 30.85613))
polyGON((52.94294 30.920592,43.61965 35.545578,55.85037 34.862696,59.12524 36.621547,47.68664 39.877048,35.69973 36.198265,37.30512 29.196711,31.09762 28.46292,52.94294 30.920592))
polyGON((52.94294 30.920592,52.45594 37.266299,59.30560 29.196711,64.12177 33.290489,58.81733 36.554277,52.94294 30.920592))

所有示例的圆心和半径:

O: (45,30)
r: 25

解决方法

对于初学者,我们可以将您的数据作为点云(多边形顶点)p[i]以及由中心p0和半径r定义的一些圆来处理。如果您的点云完全在圆内,则可以忽略半径。

我们可以利用atan2,但是为了避免交叉和扇区选择问题,我们不会像通常那样对标准笛卡尔BBOX计算扩大最小/最大范围:

  1. 计算每个点的atan2角度,并记住其在数组a[]

  2. a[]

    进行排序
  3. a[]

    中找到相应角度之间的最大距离

    别忘了角度差可以是|Pi|的顶部,因此,如果更多,则需要+/- 2*PI。还要将a[]处理为循环缓冲区。

这是我简单的C ++ / VCL尝试:

//---------------------------------------------------------------------------
float p0[]={52.87404,30.856130,42.55699,28.46292,41.54373,24.319989,53.57623,21.300564,62.94891,49.39652,27.550071,52.87404,30.85613,};
float p1[]={52.94294,30.920592,43.61965,35.545578,55.85037,34.862696,59.12524,36.621547,47.68664,39.877048,35.69973,36.198265,37.30512,29.196711,31.09762,52.94294,};
float p2[]={52.94294,52.45594,37.266299,59.30560,64.12177,33.290489,58.81733,36.554277,};
float x0=45.0,y0=30.0,R=25.0;
//---------------------------------------------------------------------------
template <class T> void sort_asc_bubble(T *a,int n)
    {
    int i,e; T a0,a1;
    for (e=1;e;n--)                                     // loop until no swap occurs
     for (e=0,a0=a[0],a1=a[1],i=1;i<n;a0=a1,i++,a1=a[i])// proces unsorted part of array
      if (a0>a1)                                        // condition if swap needed
      { a[i-1]=a1; a[i]=a0; a1=a0; e=1; }               // swap and allow to process array again
    }
//---------------------------------------------------------------------------
void get_sector(float x0,float y0,float r,float *p,int n,float &a0,float &a1)
    {
    // x0,y0    circle center
    // r        circle radius
    // p[n]     polyline vertexes
    // a0,a1    output angle range a0<=a1
    int i,j,m=n>>1;
    float x,y,*a;
    a=new float[m];
    // process points and compute angles
    for (j=0,i=0;i<n;j++)
        {
        x=p[i]-x0; i++;
        y=p[i]-y0; i++;
        a[j]=atan2(y,x);
        }
    // sort by angle
    sort_asc_bubble(a,m);
    // get max distance
    a0=a[m-1]; a1=a[0]; x=a1-a0;
    while (x<-M_PI) x+=2.0*M_PI;
    while (x>+M_PI) x-=2.0*M_PI;
    if (x<0.0) x=-x;
    for (j=1;j<m;j++)
        {
        y=a[j]-a[j-1];
        while (y<-M_PI) y+=2.0*M_PI;
        while (y>+M_PI) y-=2.0*M_PI;
        if (y<0.0) y=-y;
        if (y>x){ a0=a[j-1]; a1=a[j]; x=y; }
        }
    }
//---------------------------------------------------------------------------
void TMain::draw()
    {
    int i,n;
    float x,r,*p,a0=0.0,a1=0.0;
    float ax,ay,bx,by;
    float zoom=7.0;
    p=p0; n=sizeof(p0)/sizeof(p0[0]);
//  p=p1; n=sizeof(p1)/sizeof(p1[0]);
//  p=p2; n=sizeof(p2)/sizeof(p2[0]);

    get_sector(x0,y0,R,p,n,a0,a1);

    // clear buffer
    bmp->Canvas->Brush->Color=clBlack;
    bmp->Canvas->FillRect(TRect(0,xs,ys));

    // circle
    x=x0; y=y0; r=R;
    ax=x+R*cos(a0);
    ay=y+R*sin(a0);
    bx=x+R*cos(a1);
    by=y+R*sin(a1);
    x*=zoom; y*=zoom; r*=zoom;
    ax*=zoom; ay*=zoom;
    bx*=zoom; by*=zoom;
    bmp->Canvas->Pen->Color=clBlue;
    bmp->Canvas->Brush->Color=TColor(0x00101010);
    bmp->Canvas->Ellipse(x-r,y-r,x+r,y+r);
    bmp->Canvas->Pen->Color=clAqua;
    bmp->Canvas->Brush->Color=TColor(0x00202020);
    bmp->Canvas->Pie(x-r,y+r,ax,by);

    // PCL
    r=2.0;
    bmp->Canvas->Pen->Color=clAqua;
    bmp->Canvas->Brush->Color=clAqua;
    for (i=0;i<n;)
        {
        x=p[i]; i++;
        y=p[i]; i++;
        x*=zoom; y*=zoom;
        bmp->Canvas->Ellipse(x-r,y+r);
        }

    // render backbuffer
    Main->Canvas->Draw(0,bmp);
    }
//---------------------------------------------------------------------------

您可以忽略void TMain::draw()函数,它只是用法的一个示例,这是预览:

overview

但是,为了避免错误的结果,您有多边形(线),您有两个简单的选择:

  1. 采样点多于2点的行

    这样,角度间隙应大于点云中点之间的距离。因此,如果您对具有足够点的线进行采样,则结果将是正确的。但是,在边缘情况下,每条线选择的点数错误会导致错误的结果。另一方面,这只是添加到当前代码中的简单DDA插值。

  2. 转换为处理角度间隔,而不是角度a[]

    因此,对于每条线,请使用预定义的多边形缠绕规则(因此CW或CCW但一致)计算角度间隔<a0,a1>。而不是使用数组a[],而是对间隔列表进行了排序,在该列表中,您可以插入新的间隔,或者在重叠时与现有间隔合并。这种方法很安全,但是合并角度间隔并不是那么容易。如果输入数据是折线(如您的折线),则意味着每条下一行都从上一行端点开始,因此您可以忽略间隔列表,而只放大单个行,但是仍然需要正确处理不小的琐事。

[Edit1]使用第一种方法,更新后的功能如下:

void get_sector_pol(float x0,y0    circle center
    // r        circle radius
    // p[n]     point cloud
    // a0,k,N=10,m=(n>>1)*N;
    float ax,by,x,dx,dy,*a,_N=1.0/N;
    a=new float[m];
    // process points and compute angles
    bx=p[n-2]-x0; i++;
    by=p[n-1]-y0; i++;
    for (j=0,i=0;i<n;)
        {
        ax=bx; ay=by;
        bx=p[i]-x0; i++;
        by=p[i]-y0; i++;
        dx=_N*(bx-ax); x=ax;
        dy=_N*(by-ay); y=ay;
        for (k=0;k<N;k++,x+=dx,y+=dy,j++) a[j]=atan2(y,m);
    // get max distance
    a0=a[m-1]; a1=a[0]; x=a1-a0;
    while (x<-M_PI) x+=2.0*M_PI;
    while (x>+M_PI) x-=2.0*M_PI;
    if (x<0.0) x=-x;
    for (j=1;j<m;j++)
        {
        y=a[j]-a[j-1];
        while (y<-M_PI) y+=2.0*M_PI;
        while (y>+M_PI) y-=2.0*M_PI;
        if (y<0.0) y=-y;
        if (y>x){ a0=a[j-1]; a1=a[j]; x=y; }
        }
    }

您可以看到,它几乎相同,只是简单的 DDA 被添加到了第一行每行赢N点。在这里预览第二个多边形,该多边形由于点云方法而失败:

preview