计算几何模板

时间:2023-02-10 15:24:04

整理了一下计算几何的模板

后续会继续更新。。。。。

const double eps = 1e-8;
const double PI = acos(-1.0);
int sgn(double x){
	if(fabs(x) < eps)return 0;
	if(x < 0)return -1;
	else return 1;
}
//结构体定义 
struct Point {
	double x,y;
	Point(){}
	Point(double xx,double yy){x=xx;y=yy;}
	Point operator-(Point b){return Point(x-b.x,y-b.y);}
	Point operator+(Point b){return Point(x+b.x,y+b.y);}
	double operator^(Point b){return x*b.y-y*b.x;}
	double operator*(Point b){return x*b.x-y*b.y;}
	//绕原点逆时针转B度 
	void transXY(double B){
		double tx = x,ty = y;
		x = tx*cos(B) - ty*sin(B);
		y = tx*sin(B) + ty*cos(B);
	}
}p[maxn];
struct Line {
	Point s,e;
	Line(){}
	Line(Point ss,Point ee){s=ss;e=ee;}
	//两直线求交点 
	pair<int,Point> operator &(const Line &b)const{
		Point res = s;
		if(sgn((s-e)^(b.s-b.e)) == 0){
		if(sgn((s-b.e)^(b.s-b.e)) == 0)
		return make_pair(0,res);//重合 
		else return make_pair(1,res);//平行
		}
		double t = ((s-b.s)^(b.s-b.e))/((s-e)^(b.s-b.e));
		res.x += (e.x-s.x)*t;
		res.y += (e.y-s.y)*t;
		return make_pair(2,res);
	}
};
double dist(Point a,Point b){
	return sqrt((a-b)*(a-b));
}
bool inter(Line l1,Line l2){
	return
	max(l1.s.x,l1.e.x) >= min(l2.s.x,l2.e.x) &&
	max(l2.s.x,l2.e.x) >= min(l1.s.x,l1.e.x) &&
	max(l1.s.y,l1.e.y) >= min(l2.s.y,l2.e.y) &&
	max(l2.s.y,l2.e.y) >= min(l1.s.y,l1.e.y) &&
	sgn((l2.s-l1.e)^(l1.s-l1.e))*sgn((l2.e-l1.e)^(l1.s-l1.e)) <= 0 &&
	sgn((l1.s-l2.e)^(l2.s-l2.e))*sgn((l1.e-l2.e)^(l2.s-l2.e)) <= 0;
}
//判断直线和线段相交
bool Seg_inter_line(Line l1,Line l2) {//判断直线l1和线段l2是否相交
	return sgn((l2.s-l1.e)^(l1.s-l1.e))*sgn((l2.e-l1.e)^(l1.s-l1.e)) <= 0;
}
//点到直线距离
//返回为result,是点到直线最近的点
Point PointToLine(Point P,Line L){
	Point result;
	double t = ((P-L.s)*(L.e-L.s))/((L.e-L.s)*(L.e-L.s));
	result.x = L.s.x + (L.e.x-L.s.x)*t;
	result.y = L.s.y + (L.e.y-L.s.y)*t;
	return result;
}
//点到线段的距离
//返回点到线段最近的点
Point NearestPointToLineSeg(Point P,Line L){
	Point result;
	double t = ((P-L.s)*(L.e-L.s))/((L.e-L.s)*(L.e-L.s));
	if(t >= 0 && t <= 1){
		result.x = L.s.x + (L.e.x - L.s.x)*t;
		result.y = L.s.y + (L.e.y - L.s.y)*t;
	}
	else{
	if(dist(P,L.s) < dist(P,L.e))
	result = L.s;
	else result = L.e;
	}
	return result;
}
//计算多边形面积
//点的编号从0~n-1
double CalcArea(Point p[],int n){
	double res = 0;
	for(int i = 0;i < n;i++)
	res += (p[i]^p[(i+1)%n])/2;
	return fabs(res);
}
//*判断点在线段上
bool OnSeg(Point P,Line L){
	return
	sgn((L.s-P)^(L.e-P)) == 0 &&
	sgn((P.x - L.s.x) * (P.x - L.e.x)) <= 0 &&
	sgn((P.y - L.s.y) * (P.y - L.e.y)) <= 0;
}
//*判断点在凸多边形内
//点形成一个凸包,而且按逆时针排序(如果是顺时针把里面的<0改为>0)
//点的编号:0~n-1
//返回值:
//-1:点在凸多边形外
//0:点在凸多边形边界上
//1:点在凸多边形内
int inConvexPoly(Point a,Point p[],int n){
	for(int i = 0;i < n;i++){
		if(sgn((p[i]-a)^(p[(i+1)%n]-a)) < 0)return -1;
		else if(OnSeg(a,Line(p[i],p[(i+1)%n])))return 0;
	}
	return 1;
}
//*判断点在任意多边形内
//射线法,poly[]的顶点数要大于等于3,点的编号0~n-1
//返回值
//-1:点在凸多边形外
//0:点在凸多边形边界上
//1:点在凸多边形内
int inPoly(Point p,Point poly[],int n){
	int cnt;
	Line ray,side;
	cnt = 0;
	ray.s = p;
	ray.e.y = p.y;
	ray.e.x = -100000000000.0;//-INF,注意取值防止越界
	for(int i = 0;i < n;i++){
		side.s = poly[i];
		side.e = poly[(i+1)%n];
		if(OnSeg(p,side))return 0;
		//如果平行轴则不考虑
		if(sgn(side.s.y - side.e.y) == 0)
		continue;
		if(OnSeg(side.s,ray))
		{
		if(sgn(side.s.y - side.e.y) > 0)cnt++;
		}
		else if(OnSeg(side.e,ray))
		{
		if(sgn(side.e.y - side.s.y) > 0)cnt++;
		}
		else if(inter(ray,side))
		cnt++;
	}
	if(cnt % 2 == 1) return 1;
	else return -1;
}
//判断凸多边形
//允许共线边
//点可以是顺时针给出也可以是逆时针给出
//点的编号1~n-1
bool isconvex(Point poly[],int n){
	bool s[3];
	memset(s,false,sizeof(s));
	for(int i = 0;i < n;i++)
	{
	s[sgn( (poly[(i+1)%n]-poly[i])^(poly[(i+2)%n]-poly[i]) )+1] = true;
	if(s[0] && s[2])return false;
	}
	return true;
}
/*
* 求凸包,Graham算法
* 点的编号0~n-1
* 返回凸包结果ch[0~top-1]为凸包的编号
* poj2007 凸包点排序逆时针输出
*/
#include <stdio.h>  
#include <cstring>  
#include <cmath>  
#include <algorithm>  
using namespace std;  
  
#define sqr(a) ((a) * (a))  
#define dis(a, b) sqrt(sqr(a.x - b.x) + sqr(a.y - b.y))  
const int MAXN = 110;  
const double PI = acos(-1.0);  
struct Point {  
    int x;  
    int y;  
    Point(double a = 0, double b = 0) : x(a), y(b) {}  
    friend bool operator < (const Point &l, const Point &r) {  
        return l.y < r.y || (l.y == r.y && l.x < r.x);  
    }  
} p[MAXN], ch[MAXN];  
double mult(Point a, Point b, Point o) {  
    return (a.x - o.x) * (b.y - o.y) >= (b.x - o.x) * (a.y - o.y);  
}  
int Graham(Point p[], int n, Point res[]) {  
    int top = 1;  
    sort(p, p + n);  
    if (n == 0) return 0;  
    res[0] = p[0];  
    if (n == 1) return 0;  
    res[1] = p[1];  
    if (n == 2) return 0;  
    res[2] = p[2];  
    int i;
    for (i = 2; i < n; i++) {  
        while (top && (mult(p[i], res[top], res[top - 1])))  
            top--;  
        res[++top] = p[i];  
    }  
    int len = top;  
    res[++top] = p[n - 2];  
    for (i = n - 3; i >= 0; i--) {  
        while (top != len && (mult(p[i], res[top], res[top - 1])))  
            top--;  
        res[++top] = p[i];  
    }  
    return top;  
}  
int n;  
int main() {  
    while (scanf("%d%d", &p[n].x, &p[n].y) != EOF)  
        n++;  
    n = Graham(p, n, ch);  
    int t,i;  
    for (i = 0; i < n; i++)  
        if (ch[i].x == 0 && ch[i].y == 0) {  
            t = i;  
            break;  
        }  
    for (i = t; i < n; i++)  
        printf("(%d,%d)\n", ch[i].x, ch[i].y);  
    for (i = 0; i < t; i++)  
        printf("(%d,%d)\n", ch[i].x, ch[i].y);  
    return 0;  
}  
/*
平面最近点对 ZOJ 2107 
*/
#include<stdio.h>
#include<math.h>
#include<algorithm>
using namespace std;
const double eps=1e-8;
struct Point {
	double x,y;
}p[100005],tmp[100005];
int n;
double dist(Point a,Point b){
	return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}
int cmpxy(Point a,Point b){
	if(a.x!=b.x) return a.x<b.x;
	return a.y<b.y;
}
int cmpy(Point a,Point b){
	return a.y<b.y;
}
double Closest_Pair(int left,int right){
	double d=1e20;
	if(left==right) return d;
	if(left+1==right) return dist(p[left],p[right]);
	int mid=(left+right)>>1;
	double d1=Closest_Pair(left,mid);
	double d2=Closest_Pair(mid+1,right);
	d=min(d1,d2);
	int k=0;
	for(int i=left;i<=right;i++){
		if(fabs(p[mid].x-p[i].x)<d+eps) tmp[k++]=p[i];
	}
	sort(tmp,tmp+k,cmpy);
	for(int i=0;i<k;i++){
		for(int j=i+1;j<k&&tmp[j].y-tmp[i].y<d+eps;j++){
			d=min(d,dist(tmp[i],tmp[j]));
		}
	}
	return d;
}
int main(){
	int n;
	while(scanf("%d",&n),n){
		for(int i=0;i<n;i++) scanf("%lf%lf",&p[i].x,&p[i].y);
		sort(p,p+n,cmpxy);
		printf("%.2lf\n",Closest_Pair(0,n-1)/2);
	}
}
/*
旋转卡壳 
求平面最远点对 Beauty Contest POJ 2187
*/ 
#include <iostream>
#include <cstdlib>
#include <cstdio>
#include <cmath>
#include<algorithm>
#define N 50005
using namespace std;
const double eps=1e-5;
struct Point {
    int x,y;
    Point (){}
    Point (int xx,int yy) {x=xx;y=yy;}
    Point operator -(const Point b)const {return Point(x-b.x,y-b.y);}
    int operator ^(const Point b)const {return x*b.y-y*b.x;}
};

int n,top;
Point p[N],stack[N];
int dist(Point a,Point b){
    return (a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y);
}
int cmp(Point a,Point b) {//逆时针
    int tmp=(a-p[0])^(b-p[0]);
    if(tmp == 0 ) return dist(a,p[0]) - dist(b,p[0]) <= 0;
    else return tmp>0;
}
void Graham(){
    int k=0;
    for(int i=0;i<n;i++){
        if(p[i].y<p[k].y || (p[i].y==p[k].y && p[i].x<p[k].x)) k=i;
    }
    swap(p[0],p[k]);
    sort(p+1,p+n,cmp);
    if(n<=2){
        for(int i=0;i<n;i++) stack[i]=p[i];
        top=n;
        return;
    }
    stack[0]=p[0];  stack[1]=p[1];
    top=2;
    for(int i=2;i<n;i++){
        while(top>=2 && ((stack[top-1]-stack[top-2])^(p[i]-stack[top-2]))<=0) top--;
        stack[top++]=p[i];
    }
}
int main(){
    while(scanf("%d",&n)!=EOF){
        for(int i=0;i<n;i++) scanf("%d%d",&p[i].x,&p[i].y);
        Graham();
        for(int i=0;i<top;i++) {
            p[i]=stack[i];
        }
        p[top]=p[0];
        int ans=0;
        Point a,b;
        for(int i=0,j=1;i<top;i++){
            while(((p[i+1]-p[i])^(p[j]-p[i]))<((p[i+1]-p[i])^(p[j+1]-p[i]))) j=(j+1)%top;
              if(dist(p[i],p[j])>=ans) ans=dist(p[i],p[j]);
              if(dist(p[i+1],p[j+1])>=ans) ans=dist(p[i+1],p[j+1]);//处理两条边平行的特殊情况
        }
        printf("%d\n",ans);
    }
    return 0;
}

//旋转卡壳计算平面点集最大三角形面积
int rotating_calipers(Point p[],int n){
	int ans = 0;
	Point v;
	for(int i = 0;i < n;i++){
		int j = (i+1)%n;
		int k = (j+1)%n;
		while(j != i && k != i){
			ans = max(ans,abs((p[j]-p[i])^(p[k]-p[i])));
			while( ((p[i]-p[j])^(p[(k+1)%n]-p[k])) < 0 )
			k = (k+1)%n;
			j = (j+1)%n;
		}
	}
	return ans;
}
Point p[MAXN];
int main(){
	int n;
	while(scanf("%d",&n) == 1){
	if(n == -1)break;
	for(int i = 0;i < n;i++)list[i].input();
	Graham(n);
	for(int i = 0;i < top;i++)p[i] = list[Stack[i]];
	printf("%.2f\n",(double)rotating_calipers(p,top)/2);
	}
	return 0;
}

/*
最近点对和最远点对产生在有向平行切线上,将一条边平行过去,然后判断
求两凸包最近点对,最远点对类似 
*/ 
#include <iostream>
#include <cstdlib>
#include <cstdio>
#include <cmath>
#include<algorithm>
using namespace std;
const double eps=1e-8;
struct Point {
    double x,y;
    Point (){}
    Point (double xx,double yy) {x=xx;y=yy;}
    Point operator -(const Point b)const {return Point(x-b.x,y-b.y);}
    double operator ^(const Point b)const {return x*b.y-y*b.x;}
    double operator *(const Point &b)const{
        return x*b.x + y*b.y;
    }
};
struct Line {
    Point s,e;
    double angle;
    Line(){}
    Line(Point ss,Point ee) {s=ss;e=ee;}
};
Point p[10005],q[10005],list[10005];
int n,m,top,miny,maxy;
double dist(Point a,Point b) {
    return sqrt((b.x-a.x)*(b.x-a.x)+(b.y-a.y)*(b.y-a.y));
}
//求点到线段距离,不是点到直线距离
double pointtoseg(Point P,Point b,Point c){
    Point result;
    Line L=Line(b,c);
double t = ((P-L.s)*(L.e-L.s))/((L.e-L.s)*(L.e-L.s));
if(t >=0 && t <= 1)
{
result.x = L.s.x + (L.e.x - L.s.x)*t;
result.y = L.s.y + (L.e.y - L.s.y)*t;
}
else
{
if(dist(P,L.s) < dist(P,L.e))
result = L.s;
else result = L.e;
}
return dist(result,P);
}
double segtoseg(Point a,Point b,Point c,Point d){
    double ans1=min(pointtoseg(a,c,d),pointtoseg(b,c,d));
    double ans2=min(pointtoseg(c,a,b),pointtoseg(d,a,b));
    return min(ans1,ans2);
}

int cmp(Point a,Point b) {//逆时针
    double tmp=(a-list[0])^(b-list[0]);
    if(tmp == 0) return dist(a,list[0]) - dist(b,list[0]) <-eps;
    else return tmp>0;
}
void Graham(Point* p,Point* stack,int n){
    int k=0;
    for(int i=0;i<n;i++){
        if(p[i].y<p[k].y || (p[i].y==p[k].y && p[i].x<p[k].x)) k=i;
    }
    swap(p[0],p[k]);
    sort(p+1,p+n,cmp);
    if(n<=2){
        for(int i=0;i<n;i++) stack[i]=p[i];
        top=n;
        return;
    }
    stack[0]=p[0];  stack[1]=p[1];
    top=2;
    for(int i=2;i<n;i++){
        while(top>=2 && ((stack[top-1]-stack[top-2])^(p[i]-stack[top-2]))<eps) top--;
        stack[top++]=p[i];
    }

}
double Get_angle(Point a1,Point a2,Point b1,Point b2)
{
return (a2-a1)^(b1-b2);
}
double rotating(Point* p,int n,Point* q,int m){
    miny=0,maxy=0;
    for(int i=0;i<n;i++) {
        if(p[i].y<p[miny].y) miny=i;
    }
    for(int i=0;i<m;i++) {
        if(q[i].y>q[maxy].y) maxy=i;
    }
    double ans=dist(p[miny],q[maxy]);
    for(int i=0;i<n;i++){
        double tmp;
        while((tmp=Get_angle(p[miny],p[miny+1],q[maxy],q[maxy+1]))<-eps) maxy=(maxy+1)%m;
        if(fabs(tmp)<eps) ans=min(ans,segtoseg(p[miny],p[(miny+1)%n],q[maxy],q[(maxy+1)%m]));
        else ans=min(ans,pointtoseg(q[maxy],p[miny],p[(miny+1)%n]));
        miny=(miny+1)%n;
    }
    return ans;
}
int main() {
    while(scanf("%d%d",&n,&m)) {
        if(n==0 && m==0) break;
        //顺时针排列的点
        miny=0,maxy=0;
        for(int i=0;i<n;i++) {
            scanf("%lf%lf",&list[i].x,&list[i].y);
        }
        Graham(list,p,n);
        n=top;
        for(int i=0;i<m;i++) {
            scanf("%lf%lf",&list[i].x,&list[i].y);
        }
        Graham(list,q,m);
        m=top;
        p[n]=p[0];
        q[m]=q[0];
        printf("%.5f\n",min(rotating(p,n,q,m),rotating(q,m,p,n)));
    }
    return 0;
}

/*半平面相交(直线切割多边形)(点标号从1开始)
n^2算法 
*/
Point points[MAXN],p[MAXN],q[MAXN];
int n;
double r;
int cCnt,curCnt;
inline void getline(Point x,Point y,double &a,double &b,double &c){
    a = y.y - x.y;
    b = x.x - y.x;
    c = y.x * x.y - x.x * y.y;
}
inline void initial(){
    for(int i = 1; i <= n; ++i)p[i] = points[i];
    p[n+1] = p[1];
    p[0] = p[n];
    cCnt = n;
}
inline Point intersect(Point x,Point y,double a,double b,double c){
    double u = fabs(a * x.x + b * x.y + c);
    double v = fabs(a * y.x + b * y.y + c);
    return Point( (x.x * v + y.x * u) / (u + v) , (x.y * v + y.y * u) / (u + v) );
}
inline void cut(double a,double b ,double c){
    curCnt = 0;
    for(int i = 1; i <= cCnt; ++i){
        if(a*p[i].x + b*p[i].y + c >= EPS)q[++curCnt] = p[i];
        else {
            if(a*p[i-1].x + b*p[i-1].y + c > EPS){
                q[++curCnt] = intersect(p[i],p[i-1],a,b,c);
            }
            if(a*p[i+1].x + b*p[i+1].y + c > EPS){
                q[++curCnt] = intersect(p[i],p[i+1],a,b,c);
            }
        }
    }
    for(int i = 1; i <= curCnt; ++i)p[i] = q[i];
    p[curCnt+1] = q[1];p[0] = p[curCnt];
    cCnt = curCnt;
}
inline void solve(){
    //注意:默认点是顺时针,如果题目不是顺时针,规整化方向
    initial();
    for(int i = 1; i <= n; ++i){
        double a,b,c;
        getline(points[i],points[i+1],a,b,c);
        cut(a,b,c);
    }
    /*
    如果要向内推进r,用该部分代替上个函数
    for(int i = 1; i <= n; ++i){
        Point ta, tb, tt;
        tt.x = points[i+1].y - points[i].y;
        tt.y = points[i].x - points[i+1].x;
        double k = r / sqrt(tt.x * tt.x + tt.y * tt.y);
        tt.x = tt.x * k;
        tt.y = tt.y * k;
        ta.x = points[i].x + tt.x;
        ta.y = points[i].y + tt.y;
        tb.x = points[i+1].x + tt.x;
        tb.y = points[i+1].y + tt.y;
        double a,b,c;
        getline(ta,tb,a,b,c);
        cut(a,b,c);
    }
    */
    //多边形核的面积
    double area = 0;
    for(int i = 1; i <= curCnt; ++i)
        area += p[i].x * p[i + 1].y - p[i + 1].x * p[i].y;
    area = fabs(area / 2.0);
    //此时cCnt为最终切割得到的多边形的顶点数,p为存放顶点的数组
}
inline void GuiZhengHua(){
     //规整化方向,逆时针变顺时针,顺时针变逆时针
    for(int i = 1; i < (n+1)/2; i ++)
      swap(points[i], points[n-i]);//头文件加iostream
}
inline void init(){
     for(int i = 1; i <= n; ++i)points[i].input();
        points[n+1] = points[1];
}
/*
半平面交nlogn算法 
*/ 
#include<stdio.h>
#include<math.h>
#include<algorithm>
#define zero(a) fabs(a)<eps
using namespace std;
const double eps=1e-10;
struct Point {
    double x,y;
    Point () {};
    Point (double xx,double yy) {x=xx;y=yy;};
    double operator^(const Point b)const{
        return x*b.y-y*b.x;
    }
    Point operator-(const Point b)const{
        return Point(x-b.x,y-b.y);
    }
};
struct Line {
    Point s,e;
    double angle;
    Line(){};
    Line(Point ss,Point ee){s=ss;e=ee;}
    Point operator &(const Line b)const {
        Point res=s;
        double t=((s-b.e)^(b.s-b.e))/((s-e)^(b.s-b.e));
        res.x+=(e.x-s.x)*t;
        res.y+=(e.y-s.y)*t;
        return res;
    }
    void puts(){
        printf("%.2f %.2f %.2f %.2f  %f\n",s.x,s.y,e.x,e.y,angle);
    }
};
int n;
Line dq[20005],l[20005];
Point p[20005];
int cmp(Line a,Line b){
    if (fabs(a.angle-b.angle)<eps) return ((a.e-a.s)^(b.s-a.s))<-eps;
    else return a.angle>b.angle;
}
double xmul(Point p0,Point p1,Point p2){
    return (p1.x-p0.x)*(p2.y-p0.y)-(p2.x-p0.x)*(p1.y-p0.y);
}
void HPI(){
    sort(l,l+n,cmp);
    int ln=0;
    for(int i=1,j=0;i<n;i++) {
        if (!zero(l[i].angle-l[j].angle)) l[++j]=l[i];
        ln=j+1;
    }
    n=ln;
    int bot=0,top=1;
    dq[0]=l[0];    dq[1]=l[1];
    for(int i=2;i<n;i++){
        while(bot<top&&xmul(l[i].s,l[i].e,(dq[top]&dq[top-1]))<-eps) top--;
        while(bot<top&&xmul(l[i].s,l[i].e, (dq[bot]&dq[bot+1]))<-eps) bot++;
        dq[++top]=l[i];
    }
    while(bot<top&&xmul(dq[bot].s,dq[bot].e, (dq[top]&dq[top-1]))<-eps) top--;
    while(bot<top&&xmul(dq[top].s,dq[top].e, (dq[bot]&dq[bot+1]))<-eps) bot++;
    dq[top+1]=dq[bot];
    int cnt=0;
    for(int i=bot;i<=top;i++) p[cnt++]=dq[i]&dq[i+1];
    p[cnt]=p[0];
    if (top<=bot+1) cnt=0;
    double s=0;
    for(int i=0;i<cnt;i++) s+=p[i]^p[i+1];
    s=fabs(s*0.5);
    printf("%.1f\n",s);
}
int main() {
    while(scanf("%d",&n)!=EOF) {
        l[0]=Line(Point(0,0),Point(10000,0));   l[n].angle=atan2(l[n].e.y-l[n].s.y,l[n].e.x-l[n].s.x);
        l[0+1]=Line(Point(10000,0),Point(10000,10000)); l[0+1].angle=atan2(l[0+1].e.y-l[0+1].s.y,l[0+1].e.x-l[0+1].s.x);
        l[0+2]=Line(Point(10000,10000),Point(0,10000)); l[0+2].angle=atan2(l[0+2].e.y-l[0+2].s.y,l[0+2].e.x-l[0+2].s.x);
        l[0+3]=Line(Point(0,10000),Point(0,0)); l[0+3].angle=atan2(l[0+3].e.y-l[0+3].s.y,l[0+3].e.x-l[0+3].s.x);
        for(int i=0;i<n;i++) {
            scanf("%lf%lf%lf%lf",&l[i+4].s.x,&l[i+4].s.y,&l[i+4].e.x,&l[i+4].e.y);
            l[i+4].angle=atan2(l[i+4].e.y-l[i+4].s.y,l[i+4].e.x-l[i+4].s.x);
        }

        n=n+4;
        HPI();
    }
}

//过三点求圆心坐标
Point waixin(Point a,Point b,Point c){
	double a1 = b.x - a.x, b1 = b.y - a.y, c1 = (a1*a1 + b1*b1)/2;
	double a2 = c.x - a.x, b2 = c.y - a.y, c2 = (a2*a2 + b2*b2)/2;
	double d = a1*b2 - a2*b1;
	return Point(a.x + (c1*b2 - c2*b1)/d, a.y + (a1*c2 -a2*c1)/d);
}
//两个圆的公共部分面积
double Area_of_overlap(Point c1,double r1,Point c2,double r2){
	double d = dist(c1,c2);
	if(r1 + r2 < d + eps)return 0;
	if(d < fabs(r1 - r2) + eps){
		double r = min(r1,r2);
		return PI*r*r;
	}
	double x = (d*d + r1*r1 - r2*r2)/(2*d);
	double t1 = acos(x / r1);
	double t2 = acos((d - x)/r2);
	return r1*r1*t1 + r2*r2*t2 - d*r1*sin(t1);
}

/*
Pick公式
顶点坐标均是整点的简单多边形:面积=内部格点数目+边上格点数目/2-1
*/

/*
圆的并 
*/
/*
题目:CIRU2
题目来源:SPOJ 8119  https://www.spoj.pl/problems/CIRUT/
题目难度:中等偏难
题目内容或思路:
    圆的面积并
    题意:给n个圆(最多1000个),分别求出覆盖1层的面积、覆盖2层的面积、
        覆盖3层的面积。。。覆盖n层的面积
    方法见AC大牛blog:
    http://hi.baidu.com/aekdycoin/blog/item/c1b28e3711246b3f0b55a95e.html
做题日期:2011.3.27
*/
#include <cstdio>
#include <cstdlib>
#include <climits>
#include <iostream>
#include <algorithm>
#include <cstring>
#include <string>
#include <queue>
#include <map>
#include <vector>
#include <bitset>
#include <cmath>
#include <set>
#include <utility>
#include <ctime>
#define sqr(x) ((x)*(x))
using namespace std;

const int N = 1010;
const double eps = 1e-8;
const double pi = acos(-1.0);
double area[N];
int n;

int dcmp(double x) {
    if (x < -eps) return -1; else return x > eps;
}

struct cp {
    double x, y, r, angle;
    int d;
    cp(){}
    cp(double xx, double yy, double ang = 0, int t = 0) {
        x = xx;  y = yy;  angle = ang;  d = t;
    }
    void get() {
        scanf("%lf%lf%lf", &x, &y, &r);
        d = 1;
    }
}cir[N], tp[N * 2];

double dis(cp a, cp b) {
    return sqrt(sqr(a.x - b.x) + sqr(a.y - b.y));
}

double cross(cp p0, cp p1, cp p2) {
    return (p1.x - p0.x) * (p2.y - p0.y) - (p1.y - p0.y) * (p2.x - p0.x);
}

int CirCrossCir(cp p1, double r1, cp p2, double r2, cp &cp1, cp &cp2) {//两圆求交点数 ,并用cp1和cp2存储其交点
    double mx = p2.x - p1.x, sx = p2.x + p1.x, mx2 = mx * mx;
    double my = p2.y - p1.y, sy = p2.y + p1.y, my2 = my * my;
    double sq = mx2 + my2, d = -(sq - sqr(r1 - r2)) * (sq - sqr(r1 + r2));//sq为两圆距离的平方和
    if (d + eps < 0) return 0; if (d < eps) d = 0; else d = sqrt(d);
    double x = mx * ((r1 + r2) * (r1 - r2) + mx * sx) + sx * my2;
    double y = my * ((r1 + r2) * (r1 - r2) + my * sy) + sy * mx2;
    double dx = mx * d, dy = my * d; sq *= 2;
    cp1.x = (x - dy) / sq; cp1.y = (y + dx) / sq;
    cp2.x = (x + dy) / sq; cp2.y = (y - dx) / sq;
    if (d > eps) return 2; else return 1;
}

bool circmp(const cp& u, const cp& v) {
    return dcmp(u.r - v.r) < 0;
}

bool cmp(const cp& u, const cp& v) {
    if (dcmp(u.angle - v.angle)) return u.angle < v.angle;
    return u.d > v.d;
}

double calc(cp cir, cp cp1, cp cp2) {
    double ans = (cp2.angle - cp1.angle) * sqr(cir.r)
        - cross(cir, cp1, cp2) + cross(cp(0, 0), cp1, cp2);
    return ans / 2;
}

void CirUnion(cp cir[], int n) {
    cp cp1, cp2;
    sort(cir, cir + n, circmp);
    for (int i = 0; i < n; ++i)
        for (int j = i + 1; j < n; ++j)
            if (dcmp(dis(cir[i], cir[j]) + cir[i].r - cir[j].r) <= 0)
                cir[i].d++;//某圆被覆盖的次数
    for (int i = 0; i < n; ++i) {
        int tn = 0, cnt = 0;
        for (int j = 0; j < n; ++j) {
            if (i == j) continue;
            if (CirCrossCir(cir[i], cir[i].r, cir[j], cir[j].r,
                cp2, cp1) < 2) continue;
            cp1.angle = atan2(cp1.y - cir[i].y, cp1.x - cir[i].x);//求圆心指向交点的向量的极角,注意这里atan2(y,x)函数要先y后x
            cp2.angle = atan2(cp2.y - cir[i].y, cp2.x - cir[i].x);
            //
            cp1.d = 1;    tp[tn++] = cp1;
            cp2.d = -1;   tp[tn++] = cp2;
            if (dcmp(cp1.angle - cp2.angle) > 0) cnt++;
        }
        tp[tn++] = cp(cir[i].x - cir[i].r, cir[i].y, pi, -cnt);//第i个圆的左端点
        tp[tn++] = cp(cir[i].x - cir[i].r, cir[i].y, -pi, cnt);
        sort(tp, tp + tn, cmp);//把交点按极角排序
        int p, s = cir[i].d + tp[0].d;
        for (int j = 1; j < tn; ++j) {
            p = s;  s += tp[j].d;
            printf("%f %f %f %f\n",tp[j-1].x,tp[j].y,tp[j-1].x,tp[j].y);
            area[p] += calc(cir[i], tp[j - 1], tp[j]);//所有被覆盖了p次的圆弧的弦连起来形成的面积,等于覆盖了p次及p次以上的面积,可画图观察
        }
    }
}

void solve() {
    for (int i = 0; i < n; ++i)
        cir[i].get();
    memset(area, 0, sizeof(area));
    CirUnion(cir, n);
    for (int i = 1; i <= n; ++i) {
        area[i] -= area[i + 1];
        printf("[%d] = %.3lf\n", i, area[i]);
    }
}

int main() {
    while (scanf("%d", &n) != EOF) {
        solve();
    }
    return 0;
}