HDU 3264 - Open-air shopping malls

 

题目地址： http://acm.hdu.edu.cn/showproblem.php?pid=3264

 

计算几何。

 

思路比较简单，就是需要各种模板。。。

 

其实一开始想法就对了，不过由于二分的时候，l与r的赋值写反，一直没发现，导致后面把扫描线那题A了之后，才把这题过了。。。

 

方法很简单，先枚举以某个圆为圆心，然后再 枚举其它圆，二分寻找圆的一半面积。

 

而 覆盖面积 =  大圆弓形面积 + 小圆弓形面积  。 因为用公式求面积，需要求圆心角。

 

 求圆心角比较麻烦，我先求两圆交点，再求的两点到各自的角度。

 

两圆交点，直接从模板里找，还敲错了一个地方，囧，幸好发现了。。。

 

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<iostream>
#include<queue>
#include<map>
#include<cmath>

using namespace std;

struct Point{
	double x,y,r;
}p[30];

Point tmp;

double mx;

const double eps=1e-6;

const double sqr2 = 1.0/sqrt(2.0); 

const double pi=acos(-1.0);

int n;

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

double getangle(Point a,Point b,Point o){
	double res;
	res=atan2((a.x-o.x),(a.y-o.y))-atan2((b.x-o.x),(b.y-o.y));
	if(res<0) res+=2*pi;
	if(res>=pi) res=2*pi-res;;
	return res;
}

Point inter(Point u1,Point u2,Point v1,Point v2){
	Point ret=u1;
	double t=((u1.x-v1.x)*(v1.y-v2.y)-(u1.y-v1.y)*(v1.x-v2.x))/((u1.x-u2.x)*(v1.y-v2.y)-(u1.y-u2.y)*(v1.x-v2.x));
	ret.x+=(u2.x-u1.x)*t;
	ret.y+=(u2.y-u1.y)*t;
	return ret;
}

void inter(Point c,double r,Point l1,Point l2,Point &p1,Point &p2){
	Point p=c;
	double t;
	p.x+=l1.y-l2.y;
	p.y+=l2.x-l1.x;
	p=inter(p,c,l1,l2);
	t=sqrt(r*r-dis(p,c)*dis(p,c))/dis(l1,l2);
	p1.x=p.x+(l2.x-l1.x)*t;
	p1.y=p.y+(l2.y-l1.y)*t;
	p2.x=p.x-(l2.x-l1.x)*t;
	p2.y=p.y-(l2.y-l1.y)*t; 
}

void inter(Point c1,double r1,Point c2,double r2,Point &p1,Point &p2){
	Point u,v;
	double t;
	t=(1+(r1*r1-r2*r2)/dis(c1,c2)/dis(c1,c2))/2;
	u.x=c1.x+(c2.x-c1.x)*t;
	u.y=c1.y+(c2.y-c1.y)*t;
	v.x=u.x+c1.y-c2.y;
	v.y=u.y-c1.x+c2.x;
	inter(c1,r1,u,v,p1,p2);
} 

int main(){
	int cas,i,j;
	double l,r,mid,rdis,ag1,ag2;
	Point p1,p2,p3;
	scanf("%d",&cas);
	while(cas--){
		scanf("%d",&n);
		for(i=0;i<n;i++){
			scanf("%lf%lf%lf",&p[i].x,&p[i].y,&p[i].r);
		}
		mx=1000000;
		for(i=0;i<n;i++){
			tmp=p[i];
			tmp.r*=sqr2;
			for(j=0;j<n;j++){
				if(i==j) continue;
				rdis=dis(tmp,p[j]);
				if(rdis+p[j].r<=tmp.r) continue;
				l=rdis;
				r=rdis+p[j].r;
				while(1){
					if(r-eps<=l && l<=r+eps) break;
					mid=(l+r)/2;
					inter(tmp,mid,p[j],p[j].r,p1,p2);
					p3.x=(p1.x+p2.x)/2;
					p3.y=(p1.y+p2.y)/2;
					if(dis(p3,tmp)>=rdis){
						r=mid;
						continue;
					}
					ag1=getangle(p1,p2,tmp);
					ag2=getangle(p1,p2,p[j]);
					if(mid*mid*(ag1-sin(ag1))+p[j].r*p[j].r*(ag2-sin(ag2))>=pi*p[j].r*p[j].r){
						r=mid;
					}
					else{
						l=mid;
					}
				}
				if(l>tmp.r){
					tmp.r=l;
				}
			}
			if(tmp.r<mx){
				mx=tmp.r;
			}
		}
		printf("%.4lf\n",mx);
	}
	return 0;
}