SPOJ 8073 CIRU - The area of the union of circles（圆面积并）

Description

给出n个圆的圆心坐标和半径，求这n个圆的面积并

Input

第一行一整数n表示圆的个数，之后n行每行三个整数x[i],y[i],r[i]表示第i个圆的圆心坐标和半径(1<=n<=1000,|x[i]|,|y[i]|,r[i]<=1000)

Output

输出这n个圆的面积并，结果保留小数点后三位

Sample Input

3

0 0 1

0 0 1

100 100 1

Sample Output

6.283

Solution

首先把被包含的小圆去掉，然后把独立的圆单独计算面积，对于剩下的圆，每个圆都和某些其他圆有相交部分，把一个相交部分用一个角度区间表示（假设逆时针从0~2PI，如果一个区间[a,b]跨过了2PI，就把该区间拆成[a,2PI]和[0,b]），这些角度区间可能还有交，所以把一个圆和其他圆相交部分排序然后合并有交的区间，然后就得到了若干不相交的区间，这些区间端点都是这个圆与其他圆的交点，在圆上找到这些交点并把相邻交点连起来后，一个圆就被分成了一个凸多边形和若干圆弧，凸多边形的面积显然要被统计到答案中，而对于圆弧部分，如果两个圆相交，那么这个相交部分那块圆弧其实是在另一个圆的凸多边形部分中的，所以这种圆弧不用统计面积，而对于非相交分布的圆弧是要统计到答案中的

Code

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<vector>
#include<queue>
#include<map>
#include<set>
#include<ctime>
using namespace std;
typedef long long ll;
#define INF 0x3f3f3f3f
#define maxn 1111
#define eps 1e-8
#define PI acos(-1.0)
int sign(double x)
{
if(fabs(x)<eps)return 0;
return x>0?1:-1;
}
struct Point
{
double x,y;
Point(){};
Point(double _x,double _y){x=_x,y=_y;}
Point operator+(const Point &b)const{return Point(x+b.x,y+b.y);}
Point operator-(const Point &b)const{return Point(x-b.x,y-b.y);}
bool operator==(const Point &b)const{return sign(x-b.x)==0&&sign(y-b.y)==0;}
double norm(){return sqrt(x*x+y*y);}
};
double det(Point a,Point b)
{
return a.x*b.y-a.y*b.x;
}
double get_angle(Point a)
{
return atan2(a.y,a.x);
}
struct Reg
{
double st,ed;
Reg(){};
Reg(double _st,double _ed){st=_st,ed=_ed;}
bool operator<(const Reg&b)const
{
return sign(st-b.st)<0||sign(st-b.st)==0&&sign(ed-b.ed)<0;
}
};
struct Circle
{
double r;//半径
Point o;//圆心坐标
Circle(){};
Circle(double _r,Point _o){r=_r,o=_o;}
vector<Reg>reg;//相交部分的区间表示
Point get_point(double a){return Point(o.x+r*cos(a),o.y+r*sin(a));}
}C[maxn];
int Circle_Insert(Point o1,double r1,Point o2,double r2,Point *p)
{
double d=(o1-o2).norm();
if(sign(d)==0)
{
if(sign(r1-r2)==0)return -1;//重合
return 0;//同心圆
}
if(sign(r1+r2-d)<0)return 0;//相离
if(sign(fabs(r1-r2)-d)>0)return 0;//包含
double ang1=atan2(o2.y-o1.y,o2.x-o1.x);
double ang2=acos((r1*r1+d*d-r2*r2)/(2*r1*d));
p[0]=Point(o1.x+r1*cos(ang1+ang2),o1.y+r1*sin(ang1+ang2));
p[1]=Point(o1.x+r1*cos(ang1-ang2),o1.y+r1*sin(ang1-ang2));
if(p[0]==p[1])return 1;//相切
return 2;
}
bool mark[maxn];
double Solve(Circle *C,int n)
{
double ans=0;
memset(mark,0,sizeof(mark));
for(int i=1;i<=n;i++)
for(int j=1;j<=n;j++)
if(i!=j&&!mark[j])
{
double d=(C[i].o-C[j].o).norm();
if(sign(d+C[i].r-C[j].r)<=0)
{
mark[i]=1;
break;
}
}
for(int i=1;i<=n;i++)
if(!mark[i])
{
Point p[2];
int flag=0;
C[i].reg.clear();
for(int j=1;j<=n;j++)
if(i!=j&&!mark[j])
{
int num=Circle_Insert(C[i].o,C[i].r,C[j].o,C[j].r,p);
if(num!=2)continue;
flag=1;
double st=get_angle(p[1]-C[i].o),ed=get_angle(p[0]-C[i].o);
if(sign(st)<0)st+=2*PI;
if(sign(ed)<0)ed+=2*PI;
if(sign(st-ed)>0)
C[i].reg.push_back(Reg(st,2*PI)),C[i].reg.push_back(Reg(0,ed));
else C[i].reg.push_back(Reg(st,ed));
}
if(!flag)//与其它圆不交的圆
{
ans+=PI*C[i].r*C[i].r;
continue;
}
sort(C[i].reg.begin(),C[i].reg.end());
int cnt=1;
for(int j=1;j<C[i].reg.size();j++)
if(sign(C[i].reg[cnt-1].ed-C[i].reg[j].st)>=0)
C[i].reg[cnt-1].ed=max(C[i].reg[cnt-1].ed,C[i].reg[j].ed);
else C[i].reg[cnt++]=C[i].reg[j];
C[i].reg.push_back(C[i].reg[0]);
C[i].reg[cnt]=C[i].reg[0];
for(int j=0;j<cnt;j++)
{
p[0]=C[i].get_point(C[i].reg[j].ed);
p[1]=C[i].get_point(C[i].reg[j+1].st);
ans+=0.5*det(p[0],p[1]);
double ang=C[i].reg[j+1].st-C[i].reg[j].ed;
if(sign(ang)<0)ang+=2*PI;
ans+=0.5*C[i].r*C[i].r*(ang-sin(ang));
}
}
return ans;
}
int main()
{
int n;
scanf("%d",&n);
for(int i=1;i<=n;i++)scanf("%lf%lf%lf",&C[i].o.x,&C[i].o.y,&C[i].r);
printf("%.3f\n",Solve(C,n)+eps);
return 0;
}