HDU - 2892：area （圆与多边形交求面积）

转载

mob604756fc3573 2019-04-09 16:28:00

文章标签 #define c++ #include 凹多边形 i++ 文章分类 代码人生

pro：飞行员去轰炸一个小岛，给出炸弹落地点的位置信息，以及轰炸半径；按顺时针或者逆时针给出小岛的边界点。求被轰炸的小岛面积。

sol：即是求圆和多边形的面积交。（只会套板子的我改头换面，先理解然后手打一遍。

HDU - 2892：area （圆与多边形交求面积）_i++

首先，我们还是用三角剖分多边形，求“有向面积”，有向面积可以避免凹多边形存在的情况，或者二者相离的情况。

我们把圆心设为起点，然后对于一个三角形，求它和圆的面积交。此时有4中情况。

1，三角形在圆的内部。即OA ,OB的出度都不大于r。

2，一个点在圆外，一个点在圆内。那么交=三角形+扇形。

3，两个点在圆外，且AB与圆无交点，那么交=扇形。

4，两个点在圆外，且AB与圆有交点，那么交=三角形+扇形+三角形。

由于asin还原的区间在[-pi/2，pi/2]，所以我们在用的时候要按照角度分类讨论。

#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;
struct point{
    double x,y;
    point(){}
    point(double xx,double yy):x(xx),y(yy){}
};
struct Circle{
    point c; double r;
};
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));}
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 a[maxn];
int main()
{
    double X,Y,H,X1,Y1,t,r,ans=0; int N;
    while(~scanf("%lf%lf%lf%lf%lf%lf",&X,&Y,&H,&X1,&Y1,&r)){
       t=sqrt(H*2/10);
       X=X+X1*t; Y=Y+Y1*t; Circle C;
       C.c=point(X,Y); C.r=r;
       scanf("%d",&N);
       rep(i,1,N) scanf("%lf%lf",&a[i].x,&a[i].y);
       a[N+1]=a[1];  ans=0;
       rep(i,1,N) ans+=TriAngleCircleInsection(C,a[i],a[i+1]);
       printf("%.2lf\n",fabs(ans));
    }
    return 0;
}