您的位置:首页 > 其它

计算几何 - UVA - 12304(附2D计算几何模板)

2017-08-29 10:46 363 查看

2D几何 6合1

三角形外接圆

三角形内切圆

点到圆的切线的极角,可能有0,1,2条切线

经过给定点的圆和给定直线相切,求圆心位置,可以转化为圆与直线的交点

半径为 r 的圆与两条给定直线相切,求圆心位置,将两条直线分别延法线方向上移/下移 r 距离得到四条新直线,进而求出四个交点即为圆心位置

同时与给定的两个相离的圆外切的 半径为 r 圆的圆心,可以转化为两个半径加上 r 的大圆的交点

所有结果都要排序

问题3中极角要先计算出来再排序,不能排水平序

PS:模板是白书上的

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#include <vector>
#include <queue>
using namespace std;
typedef long long LL;

const int inf = 0x3f3f3f3f;
const double eps = 1e-6;
const double PI = acos(-1.0);

int dcmp(double x) {
if (fabs(x) < eps) return 0; else return x < 0 ? -1 : 1;
}

struct Point {
double
f695
x, y;
Point(double x = 0, double y = 0) :x(x), y(y) {}///Constructor

void debug() { printf("Point:(%.6f , %.6f)\n", x, y); }
};
typedef Point Vector;

Point read_point()
{
double x, y;
scanf("%lf%lf", &x, &y);
return Point(x, y);
}

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 * (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); }

bool operator < (const Point& a, const Point& b) {
return a.x < b.x || (a.x == b.x && a.y < b.y);
}

bool operator == (const Point& a, const Point &b) {
return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0;
}

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 Angle(Vector A, Vector B) { return acos(Dot(A, B) / Length(A) / Length(B)); }

double Cross(Vector A, Vector B) { return A.x*B.y - A.y*B.x; };
double Area2(Point A, Point B, Point C) { return Cross(B - A, C - A); }

///向量旋转,rad位弧度,逆时针
Vector Rotate(Vector A, double rad) {
return Vector(A.x*cos(rad) - A.y*sin(rad), A.x*sin(rad) + A.y*cos(rad));
}

Vector Normal(Vector A) {
double L = Length(A);
return Vector(-A.y / L, A.x / L);
}

/*得到两条直线的交点*/
Point GetLineInter(Point P, Point v, Point Q, Point w) {
Vector u = P - Q;
double t = Cross(w, u) / Cross(v, w);
return P + v*t;
}

Point GetLineProjection(Point P, Point A, Point B) {
Vector v = B - A;
return A + v*(Dot(v, P - A) / Dot(v, v));
}

double DistanceToLine(Point P, Point A, Point B) {
Vector v1 = B - A, v2 = P - A;
return fabs(Cross(v1, v2)) / Length(v1); // 如果不取绝对值,得到的是有向距离
}

struct Line {
Point p;
Vector v;

Line(Point p, Vector v) :p(p), v(v) { }

Point point(double t) {
return p + v*t;
}
Line move(double d) {
return Line(p + Normal(v)*d, v);
}
};

struct Circle {
Point c;
double r;
Circle(Point c, double r) :c(c), r(r) { }
Point point(double a) {
return Point(c.x + cos(a)*r, c.y + sin(a)*r);
}
};

Point GetLineInter(Line a, Line b) {
return GetLineInter(a.p, a.v, b.p, b.v);
}

/*输出v的极角(弧度)*/
double angle(Vector v) {
return atan2(v.y, v.x);
}

/*输出v的极角(角度),0<=ang<180*/
double lineAngleDegree(Vector v) {
double ang = angle(v)*180.0 / PI;
while (dcmp(ang) < 0) ang += 360.0;
while (dcmp(ang - 180) >= 0) ang -= 180.0;
return ang;
}

int GetLineCircleInter(Line L, Circle C, double& t1, double& t2, vector<Point>& sol) {
double a = L.v.x, b = L.p.x - C.c.x, c = L.v.y, d = L.p.y - C.c.y;
double e = a*a + c*c, f = 2 * (a*b + c*d), g = b*b + d*d - C.r*C.r;
double delta = f*f - 4 * e*g; // 判别式
if (dcmp(delta) < 0) return 0; // 相离
if (dcmp(delta) == 0) { // 相切
t1 = t2 = -f / (2 * e); sol.push_back(L.point(t1));
return 1;
}
// 相交
t1 = (-f - sqrt(delta)) / (2 * e); sol.push_back(L.point(t1));
t2 = (-f + sqrt(delta)) / (2 * e); sol.push_back(L.point(t2));
return 2;
}

int GetCircleCircleInter(Circle C1, Circle C2, vector<Point>& sol) {
double d = Length(C1.c - C2.c);
if (dcmp(d) == 0) {
if (dcmp(C1.r - C2.r) == 0) return -1; // 重合,无穷多交点
return 0;
}
if (dcmp(C1.r + C2.r - d) < 0) return 0; //相离
if (dcmp(fabs(C1.r - C2.r) - d) > 0) return 0;//包含

double a = angle(C2.c - C1.c);//c1->c2的极角
double da = acos((C1.r*C1.r + d*d - C2.r*C2.r) / (2 * C1.r*d));//差角 括号易错
Point p1 = C1.point(a - da), p2 = C1.point(a + da);

sol.push_back(p1);
if (p1 == p2) return 1;//相切
sol.push_back(p2);
return 2;//相交
}
/*点到圆的切线*/
int GetTangents(Point p, Circle C, Vector* v) {
Vector u = C.c - p;
double dist = Length(u);
if (dist < C.r) return 0;//p在圆内
else if (dcmp(dist - C.r) == 0) {// p在圆上,只有一条切线
v[0] = Rotate(u, PI / 2);
return 1;//p在圆上,只有一条切线
}
else {
double ang = asin(C.r / dist);
v[0] = Rotate(u, -ang);
v[1] = Rotate(u, +ang);
return 2;//p在圆外,有两条切线
}
}

/*两圆的切线*/
int GetTangents(Circle A, Circle B, Point* a, Point* b) {
int cnt = 0;
if (A.r < B.r) { swap(A, B);swap(a, b); }
int 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);

int rdiff = A.r - B.r;
int rsum = A.r + B.r;
if (d2 < rdiff*rdiff) return 0; //内含

double base = atan2(B.c.y - A.c.y, B.c.x - A.c.x);
if (d2 == 0 && A.r == B.r) return -1;//无限多条切线
if (d2 == rdiff*rdiff) {
a[cnt] = A.point(base);
b[cnt] = B.point(base);
cnt++;
return 1;//内切,有1条切线
}
//有外共切线
double ang = acos((A.r - B.r) / 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 (d2 == rsum*rsum) {//一条内共切线
a[cnt] = A.point(base);b[cnt] = B.point(base + PI);cnt++;
}
else if (d2>rsum*rsum) {//两条内共切线
double angg = acos((A.r + B.r) / sqrt(d2));
a[cnt] = A.point(base + angg);b[cnt] = B.point(PI + base + angg);cnt++;
a[cnt] = A.point(base - angg);b[cnt] = B.point(PI + base - angg);cnt++;
}
return cnt;
}

/*三角形外接圆*/
Circle CircumscribedCircle(Point p1, Point p2, Point p3) {
double Bx = p2.x - p1.x, By = p2.y - p1.y;
double Cx = p3.x - p1.x, Cy = p3.y - p1.y;
double D = 2 * (Bx*Cy - By*Cx);

double cx = (Cy*(Bx*Bx + By*By) - By*(Cx*Cx + Cy*Cy)) / D + p1.x;
double cy = (Bx*(Cx*Cx + Cy*Cy) - Cx*(Bx*Bx + By*By)) / D + p1.y;
Point p = Point(cx, cy);
return Circle(p, Length(p1 - p));
}

/*三角形内接圆*/
Circle InscribedCircle(Point p1, Point p2, Point p3) {
double a = Length(p2 - p3);
double b = Length(p3 - p1);
double c = Length(p1 - p2);

Point p = (p1*a + p2*b + p3*c) / (a + b + c);
return Circle(p, DistanceToLine(p, p1, p2));
}

/*Problem 4*/
vector<Point> P4(Line L, Point x, double r) {
double t1, t2;
vector<Point> sol;
GetLineCircleInter(L.move(r), Circle(x, r), t1, t2, sol);
//cout << "re:"<<re << endl;
GetLineCircleInter(L.move(-r), Circle(x, r), t1, t2, sol);
//cout << "re:" << re << endl;
return sol;
}
/*Problem 5*/
vector<Point> P5(Line a, Line b, double r) {
Line a1 = a.move(r), a2 = a.move(-r),
b1 = b.move(r), b2 = b.move(-r);
vector<Point> ans;
ans.push_back(GetLineInter(a1, b1));
ans.push_back(GetLineInter(a1, b2));
ans.push_back(GetLineInter(a2, b1));
ans.push_back(GetLineInter(a2, b2));

return ans;
}
/*Problem 6*/

vector<Point> P6(Circle a, Circle b) {
vector<Point> re;
Vector v = a.c - a.c;
double dist = Length(v);
int d = dcmp(dist - a.r - b.r);
if (d > 0) return re;
GetCircleCircleInter(a, b, re);
return re;
}

char st[100];
char q[6][100] = {
"CircumscribedCircle",
"InscribedCircle",
"TangentLineThroughPoint",
"CircleThroughAPointAndTangentToALineWithRadius",
"CircleTangentToTwoLinesWithRadius",
"CircleTangentToTwoDisjointCirclesWithRadius"
};

int main()
{
while (~scanf("%s", st)) {
if (!strcmp(st, q[0])) {
Point a = read_point(),
b = read_point(),
c = read_point();
Circle re = CircumscribedCircle(a, b, c);
printf("(%.6f,%.6f,%.6f)\n", re.c.x, re.c.y, re.r);
}
else if (!strcmp(st, q[1])) {
Point a = read_point(),
b = read_point(),
c = read_point();
Circle re = InscribedCircle(a, b, c);
printf("(%.6f,%.6f,%.6f)\n", re.c.x, re.c.y, re.r);
}
else if (!strcmp(st, q[2])) {
Point c = read_point();
double r;scanf("%lf", &r);
Point p = read_point();
vector<double> ans;
Vector v[2];
int re = GetTangents(p, Circle(c, r), v);
for (int i = 0;i<re;i++) {
ans.push_back(lineAngleDegree(v[i]));
}
sort(ans.begin(), ans.end());
printf("[");
if (re > 0)printf("%.6f", ans[0]);
if (re > 1)printf(",%.6f", ans[1]);
printf("]\n");
}
else if (!strcmp(st, q[3])) {
Point c = read_point();
Point a = read_point();
Point b = read_point();
Line L(a, b - a);
double r;scanf("%lf", &r);
vector<Point> ans;
ans = P4(L, c, r);
sort(ans.begin(), ans.end());
printf("[");
if (ans.size() > 0)printf("(%.6f,%.6f)", ans[0].x, ans[0].y);
if (ans.size() > 1)printf(",(%.6f,%.6f)", ans[1].x, ans[1].y);
printf("]\n");
}
else if (!strcmp(st, q[4])) {
Point a = read_point(), b = read_point(), c = read_point(), d = read_point();
double r;scanf("%lf", &r);
vector<Point> re = P5(Line(a, b - a), Line(c, d - c), r);
sort(re.begin(), re.end());
printf("[");
printf("(%.6f,%.6f)", re[0].x, re[0].y);
for (int i = 1;i < 4;++i) {
printf(",(%.6f,%.6f)", re[i].x, re[i].y);
}
printf("]\n");
}
else if (!strcmp(st, q[5])) {
double r1, r2, r;
Point a = read_point();
scanf("%lf", &r1);
Point b = read_point();
scanf("%lf", &r2);
scanf("%lf", &r);

vector<Point> ans = P6(Circle(a, r1 + r), Circle(b, r2 + r));
sort(ans.begin(), ans.end());
printf("[");
if (ans.size()>0)
printf("(%.6f,%.6f)", ans[0].x, ans[0].y);
if (ans.size()>1)
printf(",(%.6f,%.6f)", ans[1].x, ans[1].y);
printf("]\n");
}
}
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: