一年前写的代码,保存一下计算几何板子。
#include <bits/stdc++.h>
using namespace std;
template<typename Tp>
void Push(vector<Tp> &a, vector<Tp> b) {
for (int i = 0; i < b.size(); ++i)
a.push_back(b[i]);
return;
}
const double EPS = 1E-12, Pi = acos(-1);
inline int Sgn(double x) {
return (x < -EPS) ? (-1) : (x > EPS);
}
int Cmp(double a, double b) {
return Sgn(a - b);
}
struct Point {
double x, y;
Point(double x = 0, double y = 0) : x(x), y(y) {}
bool operator < (const Point &a) const {
return (Sgn(x - a.x)) ? (x < a.x) : (y < a.y);
}
bool operator == (const Point &a) const {
return !Sgn(x - a.x) && !Sgn(y - a.y);
}
};
typedef Point Vector;
Point Read_Point() {
Point a;
scanf("%lf%lf", &a.x, &a.y);
return a;
}
double Polar_Angle(Point a) {
return atan2(a.y, a.x);
}
double Distance(Point a, Point b) {
return sqrt(pow(a.x - b.x, 2) + pow(a.y - b.y, 2));
}
Vector operator + (Vector a, Vector b) {
return Vector(a.x + b.x, a.y + b.y);
}
Vector operator - (Point a, Point b) {
return Vector(a.x - b.x, a.y - b.y);
}
Vector operator - (Point a) {
return Vector(-a.x, -a.y);
}
Vector operator * (Vector a, double k) {
return Vector(a.x * k, a.y * k);
}
Vector operator * (double k, Vector a) {
return Vector(a.x * k, a.y * k);
}
Vector operator / (Vector a, double k) {
return Vector(a.x / k, a.y / k);
}
Vector operator / (double k, Vector a) {
return Vector(a.x / k, a.y / k);
}
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 - b.x * a.y;
}
Vector Reflect(Vector a) {
return Vector(a.x, -a.y);
}
Vector Unit(Vector a) {
return a / Length(a);
}
Vector Rotate(Vector a, Vector b) {
b = Unit(b);
return Vector(a.x * b.x - a.y * b.y, a.x * b.y + b.x * a.y);
}
Vector Normal(Vector a) {
return Vector(-a.y, a.x) / Length(a);
}
struct Angle {
Vector u, v1, v2;
Angle(Vector u, Vector v1, Vector v2) : u(u), v1(Unit(v1)), v2(Unit(v2)) {}
};
double Cosine(Angle a) {
return Dot(a.v1, a.v2);
}
struct Line {
Vector u, v;
Line(Vector u, Vector v) : u(u), v(Unit(v)) {}
};
Line Read_Line() {
Point a = Read_Point(), b = Read_Point();
return Line(a, b - a);
}
Line Rotate(Line a, Vector b) {
return Line(a.u, Rotate(a.v, b));
}
double Polar_Angle(Line a) {
if (a.v.y <= 0) a.v = -a.v;
return Polar_Angle(a.v);
}
Point Get_Line_Intersection(Line a, Line b) {
return a.u + a.v * (Cross(b.v, a.u - b.u) / Cross(a.v, b.v));
}
double Distance_To_Line(Point a, Line b) {
return fabs(Cross(b.v, a - b.u));
}
Point Get_Line_Projection(Point a, Line b) {
return b.u + b.v * Dot(b.v, a - b.u);
}
struct Segment {
Point u, v;
Segment(Point u, Point v) : u(u), v(v) {}
};
Segment Read_Segment() {
Point a = Read_Point(), b = Read_Point();
return Segment(a, b);
}
double Length(Segment a) {
return Length(a.v - a.u);
}
Point Midpoint(Segment a) {
return (a.u + a.v) / 2;
}
// double Distance_To_Segment(Point b, Segment b) {}
// bool Segment_Proper_Intersection(Segment a, Segment b) {}
// bool On_Segment(Point a, Segment b) {}
Line Angle_Bisector(Angle a) {
return Line(a.u, a.v1 + a.v2);
}
Line Perpendicular_Bisector(Segment a) {
return Line(Midpoint(a), Normal(a.v - a.u));
}
struct Circle {
Point c;
double r;
Circle(Point c, double r = 0) : c(c), r(r) {}
};
Circle Read_Circle() {
Point a = Read_Point();
double r;
scanf("%lf", &r);
return Circle(a, r);
}
Circle Circumscribed_Circle(Point a, Point b, Point c) {
Point cc = Get_Line_Intersection(Perpendicular_Bisector(Segment(a, b)), Perpendicular_Bisector(Segment(b, c)));
return Circle(cc, Distance(cc, a));
}
Circle Inscribed_Circle(Point a, Point b, Point c) {
Point cc = Get_Line_Intersection(Angle_Bisector(Angle(a, b - a, c - a)), Angle_Bisector(Angle(b, a - b, c - b)));
return Circle(cc, Distance_To_Line(cc, Line(a, b - a)));
}
vector<Point> Get_Line_Circle_Intersection(Line a, Circle b) {
vector<Point> res;
double h = Distance_To_Line(b.c, a);
if (Cmp(b.r, h) < 0) return res;
Point c = Get_Line_Projection(b.c, a);
if (!Cmp(h, b.r)) {
res.push_back(c);
return res;
}
double w = sqrt(b.r * b.r - h * h);
res.push_back(c + a.v * w);
res.push_back(c - a.v * w);
sort(res.begin(), res.end());
return res;
}
vector<Point> Get_Circle_Intersection(Circle a, Circle b) {
vector<Point> res;
if (a.r < b.r) swap(a, b);
double h = Distance(a.c, b.c);
if (Cmp(h, a.r + b.r) > 0) return res;
if (Cmp(h, a.r - b.r) < 0) return res;
if (!Cmp(h, a.r + b.r) || !Cmp(h, a.r - b.r)) {
Line x(a.c, b.c - a.c);
res.push_back(x.u + x.v * a.r);
return res;
}
double cosine = (h * h + a.r * a.r - b.r * b.r) / (2 * h * a.r);
Vector alpha(cosine, sqrt(1 - cosine * cosine)), k = Unit(b.c - a.c);
res.push_back(a.c + Rotate(k, alpha) * a.r);
res.push_back(a.c + Rotate(k, Reflect(alpha)) * a.r);
sort(res.begin(), res.end());
return res;
}
vector<double> Tangent_Line_Through_Point(Point a, Circle b) {
vector<double> res;
double h = Distance(a, b.c);
if (Cmp(h, b.r) < 0) return res;
if (!Cmp(h, b.r)) {
res.push_back(Polar_Angle(Normal(a - b.c)));
return res;
}
double w = sqrt(h * h - b.r * b.r);
Vector alpha(w, b.r);
Line x(a, b.c - a);
res.push_back(Polar_Angle(Rotate(x, alpha)) * 180 / Pi);
res.push_back(Polar_Angle(Rotate(x, Reflect(alpha))) * 180 / Pi);
sort(res.begin(), res.end());
return res;
}
vector<Point> Circle_Through_A_Point_And_Tangent_To_A_Line_With_Radius(Point a, Line b, double r) {
vector<Point> res;
Vector c = Normal(b.v);
Circle k(a, r);
Push(res, Get_Line_Circle_Intersection(Line(b.u + c * r, b.v), k));
Push(res, Get_Line_Circle_Intersection(Line(b.u - c * r, b.v), k));
sort(res.begin(), res.end());
return res;
}
vector<Point> Circle_Tangent_To_Two_Lines_With_Radius(Line a, Line b, double r) {
vector<Point> res;
Point c = Get_Line_Intersection(a, b);
Angle alpha(c, a.v, b.v);
Line x = Angle_Bisector(alpha);
double cosine = Cosine(Angle(c, x.v, b.v));
Point cc = x.u + x.v * (r / sqrt(1 - cosine * cosine));
res.push_back(cc);
res.push_back(c * 2 - cc);
alpha = Angle(c, a.v, -b.v);
x = Angle_Bisector(alpha);
cosine = Cosine(Angle(c, x.v, b.v));
cc = x.u + x.v * (r / sqrt(1 - cosine * cosine));
res.push_back(cc);
res.push_back(c * 2 - cc);
sort(res.begin(), res.end());
return res;
}
vector<Point> Circle_Tangent_To_Two_Disjoint_Circles_With_Radius(Circle a, Circle b, double r) {
return Get_Circle_Intersection(Circle(a.c, a.r + r), Circle(b.c, b.r + r));
}
char buff[100];
void Write(Point a) {
printf("(%.6f,%.6f)", a.x, a.y);
return;
}
void Write(Circle a) {
printf("(%.6f,%.6f,%.6f)", a.c.x, a.c.y, a.r);
return;
}
void Write(vector<double> a) {
putchar('[');
if (a.size()) printf("%.6f", a[0]);
for (int i = 1; i < a.size(); ++i)
printf(",%.6f", a[i]);
putchar(']');
return;
}
void Write(vector<Point> a) {
putchar('[');
if (a.size()) Write(a[0]);
for (int i = 1; i < a.size(); ++i) {
putchar(',');
Write(a[i]);
}
putchar(']');
return;
}
int main() {
while (scanf("%s", buff + 1) == 1) {
int len = strlen(buff + 1);
if (len == 19) {
Point a = Read_Point(), b = Read_Point(), c = Read_Point();
Write(Circumscribed_Circle(a, b, c));
} else if (len == 15) {
Point a = Read_Point(), b = Read_Point(), c = Read_Point();
Write(Inscribed_Circle(a, b, c));
} else if (len == 23) {
Circle a = Read_Circle();
Point b = Read_Point();
Write(Tangent_Line_Through_Point(b, a));
} else if (len == 46) {
Point a = Read_Point();
Line b = Read_Line();
double r;
scanf("%lf", &r);
Write(Circle_Through_A_Point_And_Tangent_To_A_Line_With_Radius(a, b, r));
} else if (len == 33) {
Line a = Read_Line(), b = Read_Line();
double r;
scanf("%lf", &r);
Write(Circle_Tangent_To_Two_Lines_With_Radius(a, b, r));
} else {
Circle a = Read_Circle(), b = Read_Circle();
double r;
scanf("%lf", &r);
Write(Circle_Tangent_To_Two_Disjoint_Circles_With_Radius(a, b, r));
}
putchar('
');
}
return 0;
}