【BZOJ-2618】凸多边形 计算几何 + 半平面交 + 增量法 + 三角剖分

2618: [Cqoi2006]凸多边形

Time Limit: 5 Sec  Memory Limit: 128 MB
Submit: 959  Solved: 489
[Submit][Status][Discuss]

Description

逆时针给出n个凸多边形的顶点坐标,求它们交的面积。例如n=2时,两个凸多边形如下图:
 

 

则相交部分的面积为5.233。

Input

第一行有一个整数n,表示凸多边形的个数,以下依次描述各个多边形。第i个多边形的第一行包含一个整数mi,表示多边形的边数,以下mi行每行两个整数,逆时针给出各个顶点的坐标。

Output

    输出文件仅包含一个实数,表示相交部分的面积,保留三位小数。

Sample Input

2
6
-2 0
-1 -2
1 -2
2 0
1 2
-1 2
4
0 -3
1 -1
2 2
-1 0

Sample Output

5.233

HINT

100%的数据满足:2<=n<=10,3<=mi<=50,每维坐标为[-1000,1000]内的整数

Source

Solution

裸半平面交

这里用的O(n^{2})的增量法,求完半平面交后再求多边形面积就好了,三角剖分一下就可以

一个不错的讲解

具体的做法:

•初始化时加上一个范围巨大的“框”
•每次拿一个新的半平面切割原先的凸集
•保留在新加直线左边的点,删除右边的,有向直线与多边形相交产生的新的点加入到新多边形内

Code

#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cstring>
#include<cmath>
#include<vector>
#include<cstdlib>
using namespace std;
struct Vector
{
    double x,y;
    Vector(double X=0,double Y=0) {x=X,y=Y;}
};
typedef Vector Point;
typedef vector<Point> Polygon;
Polygon polygon;
#define MAXN 20
#define MAXM 100
Point P[MAXN][MAXM];
#define eps 1e-8
#define INF 1000
const double pi= acos(-1.0);
Vector operator + (Vector A,Vector B) {return ((Vector){A.x+B.x,A.y+B.y});}
Vector operator - (Vector A,Vector B) {return ((Vector){A.x-B.x,A.y-B.y});}
Vector operator * (Vector A,double p) {return ((Vector){A.x*p,A.y*p});}
Vector operator / (Vector A,double p) {return ((Vector){A.x/p,A.y/p});}
int dcmp(double x) {if(fabs(x)<eps) return 0; else return x<0? -1:1;}
bool operator == (const Vector& a,const Vector& b) {return dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)==0;}
double Dot(Vector A,Vector B) {return A.x*B.x+A.y*B.y;}
double Len(Vector A) {return sqrt(Dot(A,A));}
double Cross(Vector A,Vector B) {return A.x*B.y-A.y*B.x;}
Point GLI(Point P,Vector v,Point Q,Vector w) {Vector u=P-Q; double t=Cross(w,u)/Cross(v,w); return P+v*t;}
double DisTL(Point P,Point A,Point B) {Vector v1=B-A,v2=P-A; return fabs(Cross(v1,v2)/Len(v1));}
bool OnSegment(Point P,Point A,Point B) {return dcmp(DisTL(P,A,B))==0&&dcmp(Dot(P-A,P-B))<0&&!(P==A)&&!(P==B);}
double PolygonArea(Polygon p)
{
    double area=0;
    int n=p.size();
    for(int i=1;i<n-1;i++)
        area+=Cross(p[i]-p[0],p[i+1]-p[0]);
    return area/2;
}
Polygon CutPolygon(Polygon poly,Point A,Point B)
{
    Polygon newpoly;
    Point C,D,ip;
    int n=poly.size(),i;
    for(i=0;i<n;i++)
    {
        C=poly[i];
        D=poly[(i+1)%n];
        if(dcmp(Cross(B-A,C-A))>=0)
            newpoly.push_back(C);
        if(dcmp(Cross(B-A,D-C))!=0)
        {
            ip=GLI(A,B-A,C,D-C);
            if(OnSegment(ip,C,D))
                newpoly.push_back(ip);
        }
    }
    return newpoly;
}
void InitPolygon(Polygon &poly,double inf)
{
    poly.clear();
    poly.push_back((Point){-inf,-inf});
    poly.push_back((Point){inf,-inf});
    poly.push_back((Point){inf,inf});
    poly.push_back((Point){-inf,inf});
}    
void Debug()
{
    for (int j=0; j<polygon.size(); j++)
        printf("(%.1lf,%.1lf)-->",polygon[j].x,polygon[j].y);
    printf("(%.1lf,%.1lf)",polygon[0].x,polygon[0].y);
    puts("");
}
int main()
{
    int N,M;
      scanf("%d",&N);
    InitPolygon(polygon,INF);
    for (int i=1; i<=N; i++)
        {            
            scanf("%d",&M);
            for (int j=1; j<=M; j++)
                scanf("%lf%lf",&P[i][j].x,&P[i][j].y);
            P[i][M+1]=P[i][1];
            for (int j=1; j<=M; j++)
                {polygon=CutPolygon(polygon,P[i][j],P[i][j+1]); /*Debug();*/}
        }
    printf("%.3lf\n",PolygonArea(polygon));
    return 0;
}

自己调个模板,p事怎么这么多QAQ

posted @ 2016-07-18 20:43  DaD3zZ  阅读(783)  评论(1编辑  收藏  举报