题意:给你n个圆,让你找出公共的正交圆。圆的正交就是两圆的圆心分别和交点连线,夹角为90度。
我们可以设出要找的圆,圆心为O,然后我们根据直角三角形的规则列出方程式,最后我们可以得知要找的圆的圆心在两圆连线的某条垂线上,假设两圆的圆心分别为a,b,过圆心O做垂线,在ab上垂足为m,可以得到am^2-bm^2=r1^2-r2^2(为定值),所以满足条件的圆的圆心在Om那条垂线上。后面就可以根据n个圆解题啦。注意细节。
#include<stdio.h>
#include<algorithm>
#include<iostream>
#include<math.h>
#define sqr(x) ((x)*(x))
using namespace std;
const int maxn = 100000+10;
const double eps = 1e-4;
const double pi = acos(-1.0);
double dcmp(double x)
{
if(fabs(x)<eps)return 0;
return x>0?1:-1;
}
double add(double a,double b)
{
// if(fabs(a+b)<eps*(fabs(a)+fabs(b)))
// return 0;
return a+b;
}
struct P
{
double x,y;
P(double x=0,double y=0):x(x),y(y)
{}
P operator +(P p)
{
return P(add(x,p.x),add(y,p.y));
}
P operator -(P p)
{
return P(add(x,-p.x),add(y,-p.y));
}
P operator *(double d)
{
return P(x*d,y*d);
}
double dot(P p)
{
return add(x*p.x,y*p.y);
}
double det(P p)
{
return add(x*p.y,-y*p.x);
}
bool operator <(const P &a)const
{
return x<a.x||(x==a.x&&y<a.y);
}
bool operator==(const P& a)const
{
return dcmp(x-a.x)==0&&dcmp(y-a.y)==0;
}
} node[maxn];
struct Circle
{
P c;
double r;
Circle() {}
Circle(P c,double r=0):c(c),r(r) {}
P point(double a)
{
return P(c.x+cos(a)*r,c.y+sin(a)*r);
}
bool operator<(const Circle &a)const
{
if(c==a.c)return r<a.r;
else if(dcmp(c.x-a.c.x)==0)return c.y<a.c.y;
return c.x<a.c.x;
}
} circle[maxn];
struct Line
{
P s,e;
Line() {}
Line(P s,P e):s(s),e(e) {}
P point(double a)
{
return s+(e-s)*a;
}
} line[maxn];
P Rotate(P a,double rad)
{
return P(a.x*cos(rad)-a.y*sin(rad),a.x*sin(rad)+a.y*cos(rad));
}
P intersection(P p1,P p2,P q1,P q2)
{
return p1+(p2-p1)*((q2-q1).det(q1-p1)/(q2-q1).det(p2-p1));
}
double dist(P a,P b)
{
return sqrt(sqr(a.x-b.x)+sqr(a.y-b.y));
}
int main()
{
int n;
while(~scanf("%d",&n))
{
int bj=1;
for(int i=0; i<n; i++)
{
scanf("%lf%lf%lf",&circle[i].c.x,&circle[i].c.y,&circle[i].r);
}
if(n==1)
{
printf("-2\n");
continue;
}
sort(circle,circle+n);
int cnt=0;
int yes=0;
for(int i=0; i<n-1; i++)
{
if(circle[i].c==circle[i+1].c&& dcmp(circle[i].r-circle[i+1].r)==0)
continue;
else if(circle[i].c==circle[i+1].c&&dcmp(circle[i].r-circle[i+1].r)!=0)
{
yes=1;
break;
}
else
{
Line AB(circle[i].c,circle[i+1].c);
double dis = dist(circle[i].c,circle[i+1].c);
// printf("%lf %lf %lf\n",circle[i].c.x,circle[i].c.y,dis);
double x=(dis*dis+circle[i].r*circle[i].r-circle[i+1].r*circle[i+1].r)/(2*dis*dis);
P p=AB.point(x);
P temp = circle[i].c-circle[i+1].c;
P q = P(-temp.y*100,temp.x*100)+p;
line[cnt++]=Line(p,q);
// printf("%lf %lf %lf %lf\n",p.x,p.y,q.x,q.y);
}
}
int cnt1=0;
int ok=1;
if(yes)
{
printf("-1\n");
continue;
}
else
{
int zh=0;
for(int i=0; i<cnt-1; i++)
{
if(dcmp((line[i].s-line[i].e).det(line[i+1].s-line[i+1].e))==0)
{
if(dcmp((line[i].s-line[i+1].s).det(line[i].s-line[i+1].e))==0&&dcmp((line[i].e-line[i+1].s).det(line[i].e-line[i+1].e))==0)
continue;
zh=1;
break;
}
node[cnt1++]=intersection(line[i].s,line[i].e,line[i+1].s,line[i+1].e);
}
if(zh)
{
printf("-1\n");
continue;
}
if(cnt1==0)
{
printf("-2\n");
continue;
}
else
{
int bj0=1;
for(int i=0; i<cnt1-1; i++)
{
if(node[i]==node[i+1])
{
continue;
}
bj0=0;
break;
}
int bj1=1;
if(!bj0)
printf("-1\n");
else
{
// printf("%lf %lf\n",node[0].x,node[0].y);
for(int i=0; i<n; i++)
if(dcmp(dist(circle[i].c,node[0])-circle[i].r)<=0)
{
bj1=0;
break;
}
if(!bj1)
printf("-1\n");
else printf("%.10f %.10f %.10f\n",node[0].x,node[0].y,sqrt(sqr(dist(circle[0].c,node[0]))-circle[0].r*circle[0].r));
}
}
}
}
return 0;
}
/*
2
0 0 5
10 0 5
*/