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

时间:2022-01-09 18:22:20

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=;
const double eps=1e-;
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 ; if(x<) return -; return ;
}
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))==) return ; //,三点一线,不构成三角形
if(dcmp(DOA-C.r)<&&dcmp(DOB-C.r)<) 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*(-x/DAB)*/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*(-y/DAB)*/r/DOB)*r*r*0.5+TS*y/DAB;
}
else if(fabs(det(OA,OB))>=r*DAB||dot(AB,AC)<=||dot(BA,BC)<=)//
{
if(dot(OA,OB)<){
if(det(OA,OB)<) 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*(-x/DAB)*/r/DOA)+asin(TS*(-y/DAB)*/r/DOB))*r*r*0.5 + TS*((x+y)/DAB-);
}
}
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
}
int tail,head;
void halfplaneintersect(int N)
{
sort(t+,t+N+,cmp);
int tot=;
rep(i,,N-) {
if(fabs(getangle(t[i])-getangle(t[i+]))>eps)
t[++tot]=t[i];
}
t[++tot]=t[N]; head=tail=;
rep(i,,tot){
while(tail>head+&&onright(t[i],q[tail],q[tail-])) tail--;
while(tail>head+&&onright(t[i],q[head+],q[head+])) head++;
q[++tail]=t[i];
}
while(tail>head+&&onright(t[head+],q[tail],q[tail-])) tail--;
}
point a[maxn],A,B;
int main()
{
int N,T,Ca=; Circle C; double x1,y1,x2,y2;
scanf("%d",&T);
while(T--){
scanf("%lf%d",&C.r,&N); C.c.x=, C.c.y=;
rep(i,,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(-,-); t[N].b=point(,-);
N++; t[N].a=point(,-); t[N].b=point(,);
N++; t[N].a=point(,); t[N].b=point(-,);
N++; t[N].a=point(-,); t[N].b=point(-,-);
point cr;
scanf("%lf%lf",&cr.x,&cr.y);
rep(i,,N){ //保证草莓在刀的左边。
point p=t[i].b-t[i].a;
point k=cr-t[i].a;
if(det(p,k)<) swap(t[i].a,t[i].b);
t[i].p=t[i].b-t[i].a;
}
halfplaneintersect(N);
double ans=,sum=pi*C.r*C.r;
q[tail+]=q[head+]; q[tail+]=q[head+];
rep(i,head+,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+]),llintersect(q[i+],q[i+]));
}
printf("Case %d: %.5lf%%\n",++Ca,ans*/sum);
}
return ;
}