HDU 3982 半平面交+圆与多边形面积交

Harry Potter and J.K.Rowling

Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 65536/65536 K (Java/Others)
Total Submission(s): 685    Accepted Submission(s): 210

Problem Description

In July 31st, last month, the author of the famous novel series J.K.Rowling celebrated her 46th birthday. Many friends gave their best wishes. They gathered together and shared one large beautiful cake.
Rowling had lots of friends, and she had a knife to cut the cake into many pieces. On the cake there was a cherry. After several cuts, the piece with the cherry was left for Rowling. Before she enjoyed it, she wondered how large this piece was, i.e., she wondered how much percentage of the cake the piece with the only cherry has.

 

Input

First line has an integer T, indicating the number of test cases.
T test cases follow. The first line of each test case has one number r (1 <= r <= 10000) and one integer n (0 <= n <= 2000), indicating the radius of the cake and the number of knife cuts. n lines follow. The i-th line has four numbers, x1, y1, x2, y2. (x1, y1) and (x2, y2) are two points on the i-th cut. (-10000 <= x1, y1, x2, y2 <= 10000) The last line of each case has two number x, y, indicating the coordinate(x, y) of the cherry.

Technical Specification

1. R, x1, y2, x2, y2, x, y all have two digits below decimal points.
2. The center of the cake is always on the point (0, 0).
3. The cherry was always on the cake and would not be on the knife cuts.

 

Output

For each test case, in one line output the case number and the percentage the piece with the cherry has of whole original cake, rounded to five fractional digits.

 

Sample Input

  
  

   
   1
1.00 2
-1.00 0.00 1.00 0.00
0.00 -1.00 0.00 1.00
0.50 0.50
  
  

 

Sample Output

  
  

   
   Case 1: 25.00000%
  
  

 题目意思：有一块半径为r的圆形蛋糕，其中心在原点，有一个人在某一个点（x，y），现把蛋糕按两点所在直线切蛋糕，求切了n次后，这个人所在蛋糕的面积占总面积的百分比。
很容易想到计算每一条直线与圆的交点，半平面交求出多边形，然后求多边形与圆的面积交，
HDU坑爹，eps=1e-8就wa，1e-5就ac，蛋疼。
代码：

/* ***********************************************
Author :_rabbit
Created Time :2014/5/4 15:03:55
File Name :20.cpp
************************************************ */
#pragma comment(linker, "/STACK:102400000,102400000")
#include <stdio.h>
#include <iostream>
#include <algorithm>
#include <sstream>
#include <stdlib.h>
#include <string.h>
#include <limits.h>
#include <string>
#include <time.h>
#include <math.h>
#include <queue>
#include <stack>
#include <set>
#include <map>
using namespace std;
#define INF 0x3f3f3f3f
#define eps 1e-5
#define pi acos(-1.0)
typedef long long ll;
int dcmp(double x){
    if(fabs(x)<eps)return 0;
    return x>0?1:-1;
}
struct Point{
    double x,y;
    Point(double _x=0,double _y=0){
        x=_x;y=_y;
    }
};
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,const double &p){
    return Point(a.x*p,a.y*p);
}
Point operator / (const Point &a,const 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 vecunit(Point a){
    return a/Length(a);
}
Point Normal(Point a){
    return Point(-a.y,a.x)/Length(a);
}
Point Rotate(Point a,double rad){
    return Point(a.x*cos(rad)-a.y*sin(rad),a.x*sin(rad)+a.y*cos(rad));
}
double Area2(Point a,Point b,Point c){
    return Length(Cross(b-a,c-a));
}
bool OnSegment(Point p,Point a1,Point a2){  
    return dcmp(Cross(a1-p,a2-p))==0&&dcmp(Dot(a1-p,a2-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);
    }
    bool operator < (const Line &L) const {
        return ang<L.ang;
    }
    Point point(double d){
        return p+(v*d);
    }
};
bool OnLeft(const Line &L,const Point &p){
    return Cross(L.v,p-L.p)>=0;
}
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;
}
Point GetLineIntersection(Line a,Line b){
    return GetLineIntersection(a.p,a.v,b.p,b.v);
}
vector<Point> HPI(vector<Line> L){
    int n=L.size();
    sort(L.begin(),L.end());//将所有半平面按照极角排序。
    int first,last;
    vector<Point> p(n);
    vector<Line> q(n);
    vector<Point> ans;
    q[first=last=0]=L[0];
    for(int i=1;i<n;i++){
        while(first<last&&!OnLeft(L[i],p[last-1]))last--;//删除顶部的半平面
        while(first<last&&!OnLeft(L[i],p[first]))first++;//删除底部的半平面
        q[++last]=L[i];//将当前的半平面假如双端队列顶部。
        if(fabs(Cross(q[last].v,q[last-1].v))<eps){//对于极角相同的，选择性保留一个。
            last--;
            if(OnLeft(q[last],L[i].p))q[last]=L[i];
        }
        if(first<last)p[last-1]=GetLineIntersection(q[last-1],q[last]);//计算队列顶部半平面交点。
    }
    while(first<last&&!OnLeft(q[first],p[last-1]))last--;//删除队列顶部的无用半平面。
    if(last-first<=1)return ans;//半平面退化
    p[last]=GetLineIntersection(q[last],q[first]);//计算队列顶部与首部的交点。
    for(int i=first;i<=last;i++)ans.push_back(p[i]);//将队列中的点复制。
    return ans;
}
double PolyArea(vector<Point> p){
    int n=p.size();
    double ans=0;
    for(int i=1;i<n-1;i++)
        ans+=Cross(p[i]-p[0],p[i+1]-p[0]);
    return fabs(ans)/2;
}
struct Circle    
{    
    Point c;    
    double r;    
    Circle(){}    
    Circle(Point c, double r):c(c), r(r){}    
    Point point(double a) //根据圆心角求点坐标    
    {    
        return Point(c.x+cos(a)*r, c.y+sin(a)*r);    
    }    
};   
  
bool InCircle(Point x,Circle c){  
    return dcmp(c.r-Length(c.c-x))>=0;  
}  
bool OnCircle(Point x,Circle c){  
    return dcmp(c.r-Length(c.c-x))==0;  
}  
int getSegCircleIntersection(Line L,Circle C,Point *sol){  
    Point nor=Normal(L.v);  
    Line p1=Line(C.c,nor);  
    Point ip=GetLineIntersection(p1,L);  
    double dis=Length(ip-C.c);  
    if(dcmp(dis-C.r)>0)return 0;  
    Point dxy=vecunit(L.v)*sqrt(C.r*C.r-dis*dis);  
    int ret=0;  
    sol[ret]=ip+dxy;  
    if(OnSegment(sol[ret],L.p,L.point(1)))ret++;  
    sol[ret]=ip-dxy;  
    if(OnSegment(sol[ret],L.p,L.point(1)))ret++;  
    return ret;  
}  
double SegCircleArea(Circle C,Point a,Point b){  
    double a1=angle(a-C.c);  
    double a2=angle(b-C.c);  
    double da=fabs(a1-a2);  
    if(da>pi)da=pi*2-da;  
    return dcmp(Cross(b-C.c,a-C.c))*da*C.r*C.r/2.0;  
}  
double PolyCircleArea(Circle C,Point *p,int n){  
    double ret=0;  
    Point sol[2];  
    p[n]=p[0];  
    for(int i=0;i<n;i++){  
        double t1,t2;  
        int cnt=getSegCircleIntersection(Line(p[i],p[i+1]-p[i]),C,sol);  
        if(cnt==0){  
            if(!InCircle(p[i],C)||!InCircle(p[i+1],C))ret+=SegCircleArea(C,p[i],p[i+1]);  
            else ret+=Cross(p[i+1]-C.c,p[i]-C.c)/2;  
        }  
        if(cnt==1){  
            if(InCircle(p[i],C)&&!InCircle(p[i+1],C))ret+=Cross(sol[0]-C.c,p[i]-C.c)/2,ret+=SegCircleArea(C,sol[0],p[i+1]);  
            else ret+=SegCircleArea(C,p[i],sol[0]),ret+=Cross(p[i+1]-C.c,sol[0]-C.c)/2;  
        }  
        if(cnt==2){  
            if((p[i]<p[i+1])^(sol[0]<sol[1]))swap(sol[0],sol[1]);  
            ret+=SegCircleArea(C,p[i],sol[0]);  
            ret+=Cross(sol[1]-C.c,sol[0]-C.c)/2;  
            ret+=SegCircleArea(C,sol[1],p[i+1]);  
        }  
    }  
    return fabs(ret);  
}  
pair<Point,Point> pp[2200];
Point p[100010];
int main()
{
     //freopen("data.in","r",stdin);
     //freopen("data.out","w",stdout);
     int T,n;
     cin>>T;
     double R;
     Point dd;
     for(int t=1;t<=T;t++){
         scanf("%lf%d",&R,&n);
         Point a,b;
         vector<Line> L;
         Line s;
         a=Point(-20000,-20000);b=Point(20000,-20000);s=Line(a,b-a);L.push_back(s);
         a=Point(20000,-20000);b=Point(20000,20000);s=Line(a,b-a);L.push_back(s);
         a=Point(20000,20000);b=Point(-20000,20000);s=Line(a,b-a);L.push_back(s);
         a=Point(-20000,10000);b=Point(-20000,-20000);s=Line(a,b-a);L.push_back(s);
         for(int i=0;i<n;i++)
             scanf("%lf%lf%lf%lf",&pp[i].first.x,&pp[i].first.y,&pp[i].second.x,&pp[i].second.y);
         scanf("%lf%lf",&dd.x,&dd.y);
         for(int i=0;i<n;i++){
             a=pp[i].first;
             b=pp[i].second;
             if(OnLeft(Line(a,b-a),dd))L.push_back(Line(a,b-a));
             else L.push_back(Line(b,a-b));
         }
         vector<Point> ff=HPI(L);
         n=ff.size();
         for(int i=0;i<n;i++)p[i]=ff[i];
         Circle C(Point(0,0),R);
         double ret=R*R*pi;
         double ans=PolyCircleArea(C,p,n);
         ans=100*ans/ret;
         printf("Case %d: %.5lf%%\n",t,ans);
     }
     return 0;
}