2015 湖南省赛计算几何--多边形的公共部分

来源:互联网 发布:一起走软件作弊 编辑:程序博客网 时间:2024/05/16 07:35

给定两个简单多边形,你的任务是判断二者是否有面积非空的公共部分。如下图,(a)中的两个矩形只有一条公共线段,没有公共面积。 


blob.png

                                                                        (a)                                                                   (b) 

在本题中,简单多边形是指不自交(也不会接触自身)、不含重复顶点并且相邻边不共线的多边形。 

注意:本题并不复杂,但有很多看上去正确的算法实际上暗藏缺陷,请仔细考虑各种情况。 

                


输入描述

输入包含不超过100组数据。每组数据包含两行,每个多边形占一行。多边形的格式是:第一个整数n表示顶点的个数 (3<=n<=100),接下来是n对整数(x,y) (-1000<=x,y<=1000),即多边形的各个顶点,按照逆时针顺序排列。 


输出描述

对于每组数据,如果有非空的公共部分,输出"Yes",否则输出"No"。 



输入样例

4 0 0 2 0 2 2 0 24 2 0 4 0 4 2 2 24 0 0 2 0 2 2 0 24 1 0 3 0 3 2 1 2

输出样例

Case 1: NoCase 2: Yes

听人说,湖南省赛计算几何好难,恰好有比赛重现,去写了一发,一开始拿有向线段暴力去切割多边形,结果wa了,后来发现好像又bug,于是想了。另外一种方法,使用多边形面积并与两个多边形面积比较,假如面积并比两个多边形面积小,说明两个多边形有面积交,否则没有。。顺便学了一下多边形面积并的写法,第一次写写了好久。。

贴代码:

#include<bits/stdc++.h>using namespace std;#define eps 1e-8#define pi acos(-1.0)int dcmp(double x){    if(fabs(x)<eps)return 0;    return x>0?1:-1;}struct Point{    double x,y;    Point(){}    Point(double _x,double _y):x(_x),y(_y){}};Point operator + (Point a,Point b){return Point(a.x+b.x,a.y+b.y);}Point operator - (Point a,Point b){return Point(a.x-b.x,a.y-b.y);}Point operator * (Point a,double p){return Point(a.x*p,a.y*p);}Point operator / (Point a,double p){return Point(a.x/p,a.y/p);}bool operator < (const Point &a,const Point &b){return a.x<b.x||(a.x==b.x&&a.y<b.y);}bool operator ==(const Point &a,const Point &b){return dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)==0;}double Dot(Point a,Point b){return a.x*b.x+a.y*b.y;}double Length(Point a){return sqrt(Dot(a,a));}double Angle(Point a,Point b){return acos(Dot(a,b)/Length(a)/Length(b));}double angle(Point a){return atan2(a.y,a.x);}double Cross(Point a,Point b){return a.x*b.y-a.y*b.x;}Point Normal(Point x){return Point(-x.y,x.x)/Length(x);}Point vecunit(Point x){return Point(x)/Length(x);}Point GetLineINtersection(Point p,Point v,Point q,Point w){    Point u=p-q;    double t=Cross(w,u)/Cross(v,w);    return p+v*t;}bool Onsegment(Point p,Point p1, Point p2){return dcmp(Cross(p1-p,p2-p))==0&&dcmp(Dot(p1-p,p2-p))<0;}struct Line{    Point p,v;    double ang;    Line(){}    Line(Point _p,Point _v):p(_p),v(_v){ang=atan2(v.y,v.x);}    Point point(double a){return p+(v*a);}    bool operator < (const Line &L)const{return ang<L.ang;}};Point GetLineIntersection(Line a,Line b){return GetLineINtersection(a.p,a.v,b.p,b.v);}double Area(vector<Point> p){    double ans=0;    int n=p.size();    for(int i=1;i<n-1;i++)        ans+=Cross(p[i]-p[0],p[i+1]-p[0]);    return ans/2;}double seg(Point o,Point a,Point b){if(dcmp(b.x-a.x)==0)return (o.y-a.y)/(b.y-a.y);return (o.x-a.x)/(b.x-a.x);}vector<Point> pp[100];pair<double,int> s[2000000];double PolyUnion(vector<Point>*p,int n){double ret=0;for(int i=0;i<n;i++){int sz1=p[i].size();for(int ii=0;ii<sz1;ii++){int M=0;s[M++]=make_pair(0,0);s[M++]=make_pair(1,0);Point A=p[i][ii],B=p[i][(ii+1)%sz1];for(int j=0;j<n;j++)if(i!=j){int sz2=p[j].size();for(int jj=0;jj<sz2;jj++){Point C=p[j][jj],D=p[j][(jj+1)%sz2];int c1=dcmp(Cross(B-A,C-A));int c2=dcmp(Cross(B-A,D-A));if(c1==0&&c2==0){if(dcmp(Dot(B-A,D-C))>0&&i>j){s[M++]=make_pair(seg(C,A,B),1);s[M++]=make_pair(seg(D,A,B),-1);}}else{double s1=Cross(D-C,A-C);double s2=Cross(D-C,B-C);if(c1>=0&&c2<0)s[M++]=make_pair(s1/(s1-s2),1);else if(c1<0&&c2>=0)s[M++]=make_pair(s1/(s1-s2),-1);}}}std::sort(s,s+M);double pre=std::min(std::max(s[0].first,0.0),1.0),now;double sum=0;int cov=s[0].second;for(int j=1;j<M;j++){now=std::min(std::max(s[j].first,0.0),1.0);if(!cov)sum+=now-pre;cov+=s[j].second;pre=now;}ret+=Cross(A,B)*sum;}}return ret/2;}int main(){int n,m,T=0;while(cin>>n){pp[0].clear();pp[1].clear();Point p;for(int i=0;i<n;i++){cin>>p.x>>p.y;pp[0].push_back(p);}cin>>m;for(int i=0;i<m;i++){cin>>p.x>>p.y;pp[1].push_back(p);}double t1=Area(pp[0])+Area(pp[1]);double t2=PolyUnion(pp,2);printf("Case %d: ",++T);if(dcmp(t1-t2)==0)puts("No");else puts("Yes");}return 0;}


0 0
原创粉丝点击