【算法】求两个凸多边形的交集的面积

问题:
编程计算两个凸多边形的交集多边形的面积

解法分析:
显然凸多边形的交集仍然是凸多边形,本解法的关键是首先求出交集凸多边
的全部顶点,然后利用凸多边的面积算法求解。
首先明确交集多边形的顶点组成。在本问题中,容易证明交集的顶点应该出自
下面的两类点 :
〔一〕已知的两个多边形 的顶点
〔二〕两个多边形的边的交点
第一类点的求法:显然如果一个多边形的某个顶点在另一个多边形内〔包括在边上〕则
该顶点是交集的一个顶点,因此我们只需依次对两个多边形的顶点进行判断,检验其是否包含
在另一个多边形中,即可求出第一类点。本程序中的子函数 inarea( )利用待判断点与另一
个多边形的全部边按顺序组成的三角形的面积总和与此多边形的面积进行大小比较,来判断该
点的归属,若大则在凸多边形的外部,不属于第一类点;否则,是交集的顶点。
第二类点的求法:两个多边形相交可能有两条边重合〔包括部分重合〕的情况,显然重
合的部分必是交集的一条边,这条边的顶点是已知两个多边形的顶点,即应属于第一类点,在
第一类点的求法中显然已包括了这种顶点,因此我们在求两个多边形的边的交点时,
可以不计算两条边重合时的交点,而只计算两条边相交〔不平行〕时的交点。在本程序的算
法中为了求两条边的交点,我们采用线段的参数方程表示法,即
x=a*t+b
y=c*t+d 其中0<=t<=1 ;
设线段的两个端点坐标分别是〔x0,y0〕和〔x1,y1〕,则线段的参数方程可写成:
x=〔x1–x0〕*t+x0
y=〔y1–y0〕*t+y0 其中 0<=t<=1;
据此可以求两条边的交点,详细的计算公式请见后面的源程序。
如果某个多边形的顶点在另一个多边形的边上,由上面的算法可知它将被分别计入两类
点当中,因此我们最后求得的交集多边形的顶点可能有部分是重复的,但由多边形的面积算法
可知这不会影响面积的大小,因此在确定两类点时,不需要做精细的划分。
在求三角形的面积时,用到了下面的公式:
┃ 1 x0 y0 ┃
S3 = (1/2)* ┃ 1 x1 y1 ┃
┃ 1 x2 y2 ┃
其中〔x0,y0〕〔x1,y1〕〔x2,y2〕为三角形的三个顶点坐标;并且 S3为三角形的有
向面积,之所以用有向面积是因为在确定交集多边形顶点的顺序时,要利用三角形的有向面积
进行方向判断。
源程序:
//VC调试通过,
#include<stdio.h> 
#include <stdlib.h>
#define ABS(x)  ((x)<0?-(x):(x))
 struct point
{
float x;
float y;
};
bool crosspt(point p,point p1,point p2,point q1,point q2) {
//求线段(p1,p2)和(q1,q2)的交点
float tmp,tp,tq;
tmp=(p2.x-p1.x)(q2.y-q1.y)-(q2.x-q1.x)(p2.y-p1.y);
if(tmp==0) return false; //两线段平行或相叠,此时不计算交点
tp=((p1.y-q1.y)(q2.x-q1.x)-(p1.x-q1.x)(q2.y-q1.y))/tmp;
if(tp>=0&&tp<=1) //交点在线段(p1,p2)中
{
tq=-((q1.y-p1.y)(p2.x-p1.x)-(q1.x-p1.x)(p2.y-p1.y))/tmp;
if(tq>=0&&tq<=1) //交点在线段(q1,q2)中
{ //若0<=tp,tq<=1,则两线段有交点,计算交点值
  p->x=(p2.x-p1.x)tp+p1.x;
p->y=(p2.y-p1.y)*tp+p1.y;
return true;
}
}
return false; //不相交
}
float S3(point *p1,point *p2,point *p3) {
//三角形的有向面积
return (p1->xp2->y+p3->xp1->y+p2->xp3->y -p3->xp2->y-p2->xp1->y-p1->xp3->y)/2;
}
bool pointorder(point p1,point p2,point p3) {
//本函数判断p1,p2,p3是否为逆时针顺序
//其中假设p1.y为三个参数中最小的一个
float tmps;
tmps=S3(&p1,&p2,&p3);
if(tmps>0)
return true;
else if(tmps<0) return false;
else { //此时三点共线或有重合点
if(p2.yp3.y)
return false;
else {
if((p1.x-p2.x)*(p2.x-p3.x)>=0)
return true;
else
return false;
}
}
}
void polysort(point *p,int n) {
//将多边形的顶点进行逆时针排序
point tmp;
int i,j;
for(i=2;i<=n;i++) {
if(p[1].y>p[i].y) {
tmp=p[1];
p[1]=p[i];
p[i]=tmp;
}
}
for(i=2;i<=n;i++) {
for(j=i+1;j<=n;j++) {
if(!pointorder(p[1],p[i],p[j])) {
tmp=p[i];
p[i]=p[j];
p[j]=tmp;
}
}
}
}
bool inarea(point p,point *pp,int n,float s) {
//判断点p是否属于n边形pp,其中多边形的面积为s
float tmps=0;
for(int i=1;i<=n;i++)
tmps+=ABS(S3(&p,pp+i,pp+i+1));
if(tmps>s)
return false;
else
return true;
}
float Sn(point *p,int n) {
//求多边形的面积
int i;
float s=0;
for(i=2;i<n;i++)
s+=ABS(S3(p+1,p+i,p+i+1));
return s;
}
void main() {
point poly1,poly2,*crosspoly;
float s1,s2; //纪录两个多边形的面积
int n1,n2;
int xb=1; //求交集的顶点时,用以记录下标
FILE *fp;
if((fp=fopen("input.dat","r"))==NULL) exit(0);
     fscanf(fp,"%d",&n1);
if((poly1=(point *)new point[n1+2])==NULL)
{
   printf("No enough mem!\n");
   exit(0);
}
for(int i=1;i<=n1;i++)
{
      fscanf(fp,"%f %f",&(poly1[i].x),&(poly1[i].y));
}
poly1[n1+1]=poly1[1];
     fscanf(fp,"%d",&n2);
if((poly2=(point *)new point[n2+2])==NULL)
{
   printf("No enough mem!\n");
   exit(0);
}
for(i=1;i<=n2;i++)
     fscanf(fp,"%f %f",&(poly2[i].x),&(poly2[i].y));
poly2[n2+1]=poly2[1];
    fclose(fp);
s1=Sn(poly1,n1);
s2=Sn(poly2,n2);
//根据我们求交集顶点的方法,可知其顶点数不超过:
i=(n1<n2?3*n1+n2:3*n2+n1);
    if((crosspoly=(point *)new point[i+1])==NULL)
{
   printf("No enough mem!\n");
   exit(0);
}
for(i=1;i<=n1;i++)                       //计算poly1中有多少个顶点
{
    if(inarea(poly1[i],poly2,n2,s2)) //被poly2所包含,即poly1
    {
        crosspoly[xb]=poly1[i];  //的顶点有多少属于交集顶点
        xb++;
    }
}
for(i=1;i<=n2;i++)                   //说明同上类似
{
    if(inarea(poly2[i],poly1,n1,s1))
    {
        crosspoly[xb]=poly2[i];
        xb++;
    }
}
int cnt;
for(i=1;i<=n1;i++)             //计算两个多边形边的交点个数
{
    cnt=0;                     //若两条边相叠(或重合),不
    for(int j=1;j<=n2;j++)     //计算其交点
    {
        if(crosspt(&crosspoly[xb],poly1[i],poly1[i+1],
                         poly2[j],poly2[j+1]))
        {
            xb++;
            cnt++;             //不包括重合相叠的情况,则
            if(cnt>1)j=n2;     //每条边和凸多边形最多有两个交点
        }
    }
}
polysort(crosspoly,xb-1);
    s1=Sn(crosspoly,xb-1);
delete poly1;
delete poly2;
delete crosspoly;
    printf("crosspoly area=%f\n",s1);
    exit(0);

}

此条目发表在程序开发分类目录,贴了, , 标签。将固定链接加入收藏夹。