二维几何模板 - 圆和球有关计算模板

圆和球有关计算模板

#include <iostream>
#include <stdio.h>
#include <string.h>
#include <algorithm>
#include <math.h>
#include <map>
#include <queue>
#include <vector>
#include <stack>
using namespace std;
const double PI = acos(-1.0);

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

/************************** 圆和球有关的计算 **************************/

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

///直线和圆的交点(函数返回的是交点的个数,参数sol存放的是交点本身)
int getLineCircleIntersection(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));
    t1 = (-f + sqrt(delta)) / (2 * e);
    sol.push_back(L.point(t2));
    return 2;
}

///两圆相交的代码(圆心分别为C1,C2,半径为r1,r2,交点为p1,p2)
int getCircleCircleIntersection(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);        ///向量C1C2的极角
    double da = acos((C1.r * C1.r + d * d - C2.r * C2.r) / (2 * C1.r * d));   ///C1C2到C1P1的角
    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;
}

///过点p到圆C的切线(v[i]是第i条切线的向量,返回切线条数)
int getTangents(Point p, Circle C, Vector* v)
{
    Vector u = C.c - p;
    double dist = Length(u);
    if(dist < C.r)
        return 0;
    else if(dcmp(dist - C.r) == 0)
    {
        ///p在圆上，只有一条切线
        v[0] = Rotate(u, PI / 2);
        return 1;
    }
    else
    {
        double ang = asin(C.r / dist);
        v[0] = Rotate(u, -ang);
        v[1] = Rotate(u, +ang);
        trturn 2;
    }
}

///两圆外公切线(假定r1>=r2,返回切线的条数。-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);
    }
    int d2 = (A.x - B.x) * (A.x - B.x) + (A.y - B.y) * (A.y - B.y);
    int rdiff = A.r - B.r;
    int rsum = A.r + B.r;
    if(d2 < rdiff * rdiff)///内含
        return 0;
    double base = atan2(B.y - A.y, B.x - A.x);
    if(d2 == 0 && A.r == B.r)///无限条切线
        return -1;
    if(d2 == rdiff * rdiff)
    {
        ///内切,1条切线
        a[cnt] = A.getPoint(base);
        b[cnt] = B.getPoint(base);
        cnt++;
        return 1;
    }
    ///有外公切线
    double ang = acos((A.r-B.r)sqrt(d2));
    a[cnt] = A.getPoint(base + ang);
    b[cnt] = B.getPoint(base + ang);
    cnt++;
    a[cnt] = A.getPoint(base - ang);
    b[cnt] = B.getPoint(base - ang);
    cnt++;
    if(d2 == rsum * rsum)///一条公切线
    {
        a[cnt] = A.getPoint(base);
        b[cnt] = B.getPoint(PI+base);
        cnt++;
    }
    else if(d2 > rsum * rsum)///两条公切线
    {
        double ang = acos((A.r + b.r) / sqrt(d2));
        a[cnt] = A.getPoint(base + ang);
        b[cnt] = B.getPoint(PI + base + ang);
        cnt++;
        a[cnt] = A.getPoint(base - ang);
        b[cnt] = B.getPoint(PI + base - ang);
        cnt++;
    }
    return cnt;
}