1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167 | #include <algorithm>
#include <cmath>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <vector>
using namespace std;
constexpr double EPS = 1e-8; // 精度系数
const double PI = acos(-1.0); // π
constexpr int N = 4;
// 点的定义
struct Point {
double x, y;
Point(double x = 0, double y = 0) : x(x), y(y) {}
bool operator<(Point A) const { return x == A.x ? y < A.y : x < A.x; }
};
// 向量的定义
using Vector = Point;
// 向量加法
Vector operator+(Vector A, Vector B) { return Vector(A.x + B.x, A.y + B.y); }
// 向量减法
Vector operator-(Vector A, Vector B) { return Vector(A.x - B.x, A.y - B.y); }
// 向量数乘
Vector operator*(Vector A, double p) { return Vector(A.x * p, A.y * p); }
// 向量数除
Vector operator/(Vector A, double p) { return Vector(A.x / p, A.y / p); }
// 与0的关系
int dcmp(double x) {
if (fabs(x) < EPS) return 0;
return x < 0 ? -1 : 1;
}
// 向量点乘
double Dot(Vector A, Vector B) { return A.x * B.x + A.y * B.y; }
// 向量长度
double Length(Vector A) { return sqrt(Dot(A, A)); }
// 向量叉乘
double Cross(Vector A, Vector B) { return A.x * B.y - A.y * B.x; }
// 点在直线上投影
Point GetLineProjection(Point P, Point A, Point B) {
Vector v = B - A;
return A + v * (Dot(v, P - A) / Dot(v, v));
}
// 圆
struct Circle {
Point c;
double r;
Circle() : c(Point(0, 0)), r(0) {}
Circle(Point c, double r = 0) : c(c), r(r) {}
// 输入极角返回点坐标
Point point(double a) { return Point(c.x + cos(a) * r, c.y + sin(a) * r); }
};
// 两圆公切线 返回切线的条数,-1表示无穷多条切线
// a[i] 和 b[i] 分别是第i条切线在圆A和圆B上的切点
int getTangents(Circle A, Circle B, Point* a, Point* b) {
int cnt = 0;
if (A.r < B.r) {
swap(A, B);
swap(a, b);
}
double d2 =
(A.c.x - B.c.x) * (A.c.x - B.c.x) + (A.c.y - B.c.y) * (A.c.y - B.c.y);
double rdiff = A.r - B.r;
double rsum = A.r + B.r;
if (dcmp(d2 - rdiff * rdiff) < 0) return 0; // 内含
double base = atan2(B.c.y - A.c.y, B.c.x - A.c.x);
if (dcmp(d2) == 0 && dcmp(A.r - B.r) == 0) return -1; // 无限多条切线
if (dcmp(d2 - rdiff * rdiff) == 0) { // 内切,一条切线
a[cnt] = A.point(base);
b[cnt] = B.point(base);
++cnt;
return 1;
}
// 有外公切线
double ang = acos(rdiff / sqrt(d2));
a[cnt] = A.point(base + ang);
b[cnt] = B.point(base + ang);
++cnt;
a[cnt] = A.point(base - ang);
b[cnt] = B.point(base - ang);
++cnt;
if (dcmp(d2 - rsum * rsum) == 0) { // 一条内公切线
a[cnt] = A.point(base);
b[cnt] = B.point(PI + base);
++cnt;
} else if (dcmp(d2 - rsum * rsum) > 0) { // 两条内公切线
double ang = acos(rsum / sqrt(d2));
a[cnt] = A.point(base + ang);
b[cnt] = B.point(PI + base + ang);
++cnt;
a[cnt] = A.point(base - ang);
b[cnt] = B.point(PI + base - ang);
++cnt;
}
return cnt;
}
// 点 O 在圆 A 外,求圆 A 的反演圆 B,R 是反演半径
Circle Inversion_C2C(Point O, double R, Circle A) {
double OA = Length(A.c - O);
double RB = 0.5 * ((1 / (OA - A.r)) - (1 / (OA + A.r))) * R * R;
double OB = OA * RB / A.r;
double Bx = O.x + (A.c.x - O.x) * OB / OA;
double By = O.y + (A.c.y - O.y) * OB / OA;
return Circle(Point(Bx, By), RB);
}
// 直线反演为过 O 点的圆 B,R 是反演半径
Circle Inversion_L2C(Point O, double R, Point A, Vector v) {
Point P = GetLineProjection(O, A, A + v);
double d = Length(O - P);
double RB = R * R / (2 * d);
Vector VB = (P - O) / d * RB;
return Circle(O + VB, RB);
}
// 返回 true 如果 A B 两点在直线同侧
bool theSameSideOfLine(Point A, Point B, Point S, Vector v) {
return dcmp(Cross(A - S, v)) * dcmp(Cross(B - S, v)) > 0;
}
int main() {
int T;
scanf("%d", &T);
while (T--) {
Circle A, B;
Point P;
scanf("%lf%lf%lf", &A.c.x, &A.c.y, &A.r);
scanf("%lf%lf%lf", &B.c.x, &B.c.y, &B.r);
scanf("%lf%lf", &P.x, &P.y);
Circle NA = Inversion_C2C(P, 10, A);
Circle NB = Inversion_C2C(P, 10, B);
Point LA[N], LB[N];
Circle ansC[N];
int q = getTangents(NA, NB, LA, LB), ans = 0;
for (int i = 0; i < q; ++i)
if (theSameSideOfLine(NA.c, NB.c, LA[i], LB[i] - LA[i])) {
if (!theSameSideOfLine(P, NA.c, LA[i], LB[i] - LA[i])) continue;
ansC[ans++] = Inversion_L2C(P, 10, LA[i], LB[i] - LA[i]);
}
printf("%d\n", ans);
for (int i = 0; i < ans; ++i) {
printf("%.8f %.8f %.8f\n", ansC[i].c.x, ansC[i].c.y, ansC[i].r);
}
}
return 0;
}
|