1. 程式人生 > >【CQOI2006】凸多邊形

【CQOI2006】凸多邊形

1713 -- 【CQOI2006】凸多邊形

Description

  逆時針給出n個凸多邊形的頂點座標,求它們交的面積。例如n=2時,兩個凸多邊形如下圖:
        img
  則相交部分的面積為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

\(\\\)

半平面交的模板。

\(atan2\)求的是奇角(值域是\((-\pi,\pi]\))。

程式碼:

#include<bits/stdc++.h>
#define ll long long
#define N 200005
#define eps 1e-12

using namespace std;
inline int Get() {int x=0,f=1;char ch=getchar();while(ch<'0'||ch>'9') {if(ch=='-') f=-1;ch=getchar();}while('0'<=ch&&ch<='9') {x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}return x*f;}

struct Point {
    double x,y;
    double angle() {return atan2(y,x);}
    void out() {cout<<"("<<x<<","<<y<<") ";}
};
typedef Point Vector;
Point operator + (const Point &a,const Point &b) {return (Point) {a.x+b.x,a.y+b.y};}
Point operator - (const Point &a,const Point &b) {return (Point) {a.x-b.x,a.y-b.y};}
Point operator * (const Point &a,double b) {return (Point) {a.x*b,a.y*b};}
Point operator / (const Point &a,double b) {return (Point) {a.x/b,a.y/b};}

double Cross(const Vector &a,const Vector &b) {return a.x*b.y-a.y*b.x;}

struct Line {
    Point s,t;
    double ang;
    void Init(Point a,Point b) {
        s=a,t=b;
        ang=(t-s).angle();
    }
    void out() {s.out(),t.out();cout<<ang;puts("");}
}l[N<<1];

int tot;
bool OnRight(const Line &a,const Point &b) {return Cross(a.t-a.s,b-a.s)<-eps;}
bool operator <(const Line &a,const Line &b) {
    double x=a.ang-b.ang;
    if(fabs(x)>eps) return x<0;
    return OnRight(a,b.t);
}

Point intersection(const Line &a,const Line &b) {
    double S=Cross(a.t-a.s,b.s-a.s),T=Cross(a.t-a.s,b.t-a.s);
    return b.s+(b.t-b.s)*(S/(S-T));
}
bool is_parallel(const Line &a,const Line &b) {
    return fabs(Cross(a.t-a.s,b.t-b.s))<eps;
}

int n,m;
Point p[N];
bool SI(Line *l,int n,int &m) {
    static Line q[N];
    static Point q2[N];
    int h=0,t=0;
    q[h=t=1]=l[1];
    for(int i=2;i<=n;i++) {
        if(fabs(l[i].ang-l[i-1].ang)<eps) continue ;
        if(h<t&&(is_parallel(q[h],q[h+1])||is_parallel(q[t],q[t-1]))) {
            return 0;
        }
        while(h<t&&OnRight(l[i],q2[t-1])) t--;
        while(h<t&&OnRight(l[i],q2[h])) h++;
        q[++t]=l[i];
        if(h<t) q2[t-1]=intersection(q[t],q[t-1]);
    }
    while(h<t&&OnRight(q[h],q2[t-1])) t--;
    while(h<t&&OnRight(q[t],q2[h])) h++;
    if(t-h<=1) return 0;
    q2[t]=intersection(q[h],q[t]);
    m=t-h+1;
    for(int i=1;i<=m;i++) p[i]=q2[i+h-1];
    return 1;
}

double area(Point *p,int n) {
    double ans=0;
    for(int i=2;i<n;i++) {
        ans+=Cross(p[i]-p[1],p[i+1]-p[1]);
    }
    return ans/2.0;
}

int main() {
    int T=Get(); 
    while(T--) {
        n=Get();
        for(int i=1;i<=n;i++) scanf("%lf%lf",&p[i].x,&p[i].y);
        p[n+1]=p[1];
        for(int i=1;i<=n;i++) l[++tot].Init(p[i],p[i+1]);
    }
    sort(l+1,l+1+tot);
    if(SI(l,tot,n)) {
        cout<<fixed<<setprecision(3)<<area(p,n);
    } else cout<<"0.000";
    return 0;
}