HDU4773 Problem of Apollonius(反演变换)
题意:
给定两个相离的圆,一个点P,求出过点P并且和两个圆都外切的所有圆。
思路:
以P为反演中心将两个圆O1,O2分别反演,得到O1',O2',则和O1,O2相切的圆反演后是一条不经过P点且与O1'、O2‘相切的直线。显然这条直线就是O1'、O2’的公切线。因为要求和O1、O2外切的圆,因此反演后的直线必定使O1'、O2'、P在圆的同一侧。筛选出符合条件的公切线再反演回去就是所求的圆了。
代码:
#include <bits/stdc++.h>
using namespace std;
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;
return 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);
}
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;
}
Point operator *(double b)const {
return Point(x*b,y*b);
}
Point Rotate(double rad){ //逆时针旋转
return Point(x*cos(rad)-y*sin(rad),x*sin(rad)+y*cos(rad));
}
double len(){//Vector
return sqrt(x*x+y*y);
}
void stdd(){//Vector
double le=len();
x/=le;
y/=le;
}
void input(){
scanf("%lf%lf",&x,&y);
}
void output(){
printf("(%.8f,%.8f) ",x,y);
}
};
struct Circle {
Point o;
double r;
Circle(){}
void input(){
o.input();
scanf("%lf",&r);
}
void output(){
printf("%.8f %.8f %.8f
",o.x,o.y,r);
}
};
double xmult(Point p0,Point p1,Point p2) { //p0p1 X p0p2
return (p1-p0)^(p2-p0);
}
double dist(Point a,Point b) {
return sqrt( (b - a)*(b - a) );
}
double getArea(Point a,Point b,Point c){//三角形面积S
return fabs(xmult(a,b,c))/2.0;
}
double distLP(Point a,Point b,Point p){//p到ab的距离
return getArea(a,b,p)*2/dist(a,b);
}
Circle CtC(Circle C,Point p,double r){//圆到圆的反演
Circle res;
double t = dist(C.o,p);
double x = r*r / (t - C.r);
double y = r*r / (t + C.r);
res.r = (x - y) / 2.0;
double s = (x + y) / 2.0;
res.o = p + (C.o - p) * (s / t);
return res;
}
Circle LtC(Point a,Point b,Point p,double r){//直线到圆的反演
double d=distLP(a,b,p);
d=r*r/d;
Circle c;
c.r=d/2;
Point v1;
if(xmult(a,b,p)>0)
v1=(a-b).Rotate(PI/2);
else
v1=(b-a).Rotate(PI/2);
v1.stdd();
c.o=p+v1*c.r;
return c;
}
Circle c[10];
Point p;
bool check(Point a,Point b,Point p){
if(sgn(xmult(a,b,p))==sgn(xmult(a,b,c[0].o))&&sgn(xmult(a,b,p))==sgn(xmult(a,b,c[1].o)))
return 1;
else return 0;
}
vector<Circle>ans;
void solve(){
ans.clear();
double r=1;
c[0]=CtC(c[0],p,r);
c[1]=CtC(c[1],p,r);
if(c[0].r>c[1].r)swap(c[0],c[1]);
Point v1=c[0].o-c[1].o;
double d=dist(c[0].o,c[1].o);
Point a,b;
double alpha=acos((c[1].r-c[0].r)/d);
v1.stdd();
a=c[1].o+v1.Rotate(-alpha)*c[1].r;
b=c[0].o+v1.Rotate(-alpha)*c[0].r;
if(check(a,b,p))ans.push_back(LtC(a,b,p,r));
a=c[1].o+v1.Rotate(alpha)*c[1].r;
b=c[0].o+v1.Rotate(alpha)*c[0].r;
if(check(a,b,p))ans.push_back(LtC(a,b,p,r));
}
int main () {
int T;
scanf("%d",&T);
while(T--){
c[0].input();c[1].input();
scanf("%lf%lf",&p.x,&p.y);
solve();
printf("%d
",ans.size());
for(int i=0;i<ans.size();i++){
ans[i].output();
}
}
}