Live2d Test Env

HDU - 3982:Harry Potter and J.K.Rowling(半平面交+圆与多边形求交)(WA ing)

pro:给定一枚蛋糕,蛋糕上某个位置有个草莓,寿星在上面切了N刀,最后寿星会吃含有草莓的那一块蛋糕,问他的蛋糕占总蛋糕的面积比。

sol:显然需要半平面交求含有蛋糕的那一块,然后有圆弧,不太方便求交。 所以我们可以直线构成的边界,求出平面交; 然后用这个多边形去和圆求交。

(百度了一下很多人都没过,好像是这题很卡精度,反正我每个地方都改过,还是WA,大概wa了4个小时了,要不以后再回来改。

当然也不排除有其他问题。

#include<bits/stdc++.h>
#define rep(i,a,b) for(int i=a;i<=b;i++)
using namespace std;
const int maxn=200010;
const double eps=1e-12;
const double pi=acos(-1.0);
struct point{
    double x,y;
    point(){}
    point(double xx,double yy):x(xx),y(yy){}
};
struct Circle{
    point c; double r;
};
struct line{
    point a,b;//起点
    point p;//起点到终点的向量
    double angle;
};
double det(point a,point b){ return a.x*b.y-a.y*b.x;}
double dot(point a,point b){ return a.x*b.x+a.y*b.y;}
point operator *(point a,double t){ return point(a.x*t,a.y*t);}
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);}
double Length(point A){return sqrt(dot(A,A));}
double getangle(point a){ return atan2(a.y,a.x);}
double getangle(line a){ return getangle(a.p);}
int dcmp(double x){
    if(fabs(x)<eps) return 0;  if(x<0) return -1; return 1;
}
double TriAngleCircleInsection(Circle C, point A, point B)
{
    point OA=A-C.c,OB=B-C.c;
    point BA=A-B, BC=C.c-B;
    point AB=B-A, AC=C.c-A;
    double DOA=Length(OA),DOB=Length(OB),DAB=Length(AB),r=C.r;
    if(dcmp(det(OA,OB))==0) return 0; //,三点一线,不构成三角形
    if(dcmp(DOA-C.r)<0&&dcmp(DOB-C.r)<0) return det(OA,OB)*0.5; //内部
    else if(DOB<r&&DOA>=r) //一内一外
    {
        double x=(dot(BA,BC)+sqrt(r*r*DAB*DAB-det(BA,BC)*det(BA,BC)))/DAB;
        double TS=det(OA,OB)*0.5;
        return asin(TS*(1-x/DAB)*2/r/DOA)*r*r*0.5+TS*x/DAB;
    }
    else if(DOB>=r&&DOA<r)// 一外一内
    {
        double y=(dot(AB,AC)+sqrt(r*r*DAB*DAB-det(AB,AC)*det(AB,AC)))/DAB;
        double TS=det(OA,OB)*0.5;
        return asin(TS*(1-y/DAB)*2/r/DOB)*r*r*0.5+TS*y/DAB;
    }
    else if(fabs(det(OA,OB))>=r*DAB||dot(AB,AC)<=0||dot(BA,BC)<=0)//
    {
        if(dot(OA,OB)<0){
            if(det(OA,OB)<0) return (-acos(-1.0)-asin(det(OA,OB)/DOA/DOB))*r*r*0.5;
            else  return ( acos(-1.0)-asin(det(OA,OB)/DOA/DOB))*r*r*0.5;
        }
        else      return asin(det(OA,OB)/DOA/DOB)*r*r*0.5; //小于90度,以为asin对应的区间是[-90度,90度]
    }
    else //弧+三角形
    {
        double x=(dot(BA,BC)+sqrt(r*r*DAB*DAB-det(BA,BC)*det(BA,BC)))/DAB;
        double y=(dot(AB,AC)+sqrt(r*r*DAB*DAB-det(AB,AC)*det(AB,AC)))/DAB;
        double TS=det(OA,OB)*0.5;
        return (asin(TS*(1-x/DAB)*2/r/DOA)+asin(TS*(1-y/DAB)*2/r/DOB))*r*r*0.5 + TS*((x+y)/DAB-1);
    }
}
point llintersect(line A,line B)
{
    point C=A.a-B.a;
    double t=det(C,B.p)/det(B.p,A.p);
    return A.a+A.p*t;
}
point s[maxn]; line t[maxn],q[maxn];
bool cmp(line a,line b){
    double A=getangle(a),B=getangle(b);
    point t=(b.a+b.p)-a.a;
    if(fabs(A-B)<eps) return det(a.p,t)>=0.0;
    return A<B;
}
bool onright(line P,line a,line b)
{
    point o=llintersect(a,b);
    point Q=o-P.a;
    return det(Q,P.p)>0; //如果同一直线上不能相互看到,则>=0
}
int tail,head;
void halfplaneintersect(int N)
{
    sort(t+1,t+N+1,cmp);
    int tot=0;
    rep(i,1,N-1) {
        if(fabs(getangle(t[i])-getangle(t[i+1]))>eps)
          t[++tot]=t[i];
    }
    t[++tot]=t[N]; head=tail=0;
    rep(i,1,tot){
        while(tail>head+1&&onright(t[i],q[tail],q[tail-1])) tail--;
        while(tail>head+1&&onright(t[i],q[head+1],q[head+2])) head++;
        q[++tail]=t[i];
    }
    while(tail>head+1&&onright(t[head+1],q[tail],q[tail-1])) tail--;
}
point a[maxn],A,B;
int main()
{
    int N,T,Ca=0; Circle C; double x1,y1,x2,y2;
    scanf("%d",&T);
    while(T--){
        scanf("%lf%d",&C.r,&N); C.c.x=0, C.c.y=0;
        rep(i,1,N) scanf("%lf%lf%lf%lf",&t[i].a.x,&t[i].a.y,&t[i].b.x,&t[i].b.y);
        N++; t[N].a=point(-10000,-10000); t[N].b=point(10000,-10000);
        N++; t[N].a=point(10000,-10000); t[N].b=point(10000,10000);
        N++; t[N].a=point(10000,10000); t[N].b=point(-10000,10000);
        N++; t[N].a=point(-10000,10000); t[N].b=point(-10000,-10000);
        point cr;
        scanf("%lf%lf",&cr.x,&cr.y);
        rep(i,1,N){ //保证草莓在刀的左边。
            point p=t[i].b-t[i].a;
            point k=cr-t[i].a;
            if(det(p,k)<0) swap(t[i].a,t[i].b);
            t[i].p=t[i].b-t[i].a;
        }
        halfplaneintersect(N);
        double ans=0,sum=pi*C.r*C.r;
        q[tail+1]=q[head+1]; q[tail+2]=q[head+2];
        rep(i,head+1,tail){
            //cout<<q[i].a.x<<" "<<q[i].a.y<<" "<<q[i].a.x+q[i].p.x<<" "<<q[i].a.y+q[i].p.y<<endl;
            ans+=TriAngleCircleInsection (C,llintersect(q[i],q[i+1]),llintersect(q[i+1],q[i+2]));
        }
        printf("Case %d: %.5lf%%\n",++Ca,ans*100/sum);
    }
    return 0;
}

 

posted @ 2019-04-09 21:45  nimphy  阅读(256)  评论(0编辑  收藏  举报