计算几何模版

此模板包含了一些基本简单的二维几何问题,

1三角形外接圆            2三角形内切圆

3过圆外某点切线的角度    4过某条直线外一点半径为r的圆

5和两条相交直线相切的半径为r的圆         6和两个相离的圆相切的圆

1.计算向量点积, 叉积, 长度, 夹角, 向量的旋转(逆时针), 向量的单位法线

2.计算两点距离, 点到直线距离,两直线交点, 点到线段距离, 点在直线的投影,将直线AB沿法线方向平移d得到的直线EF

3. 圆与直线的交点(相离,没有交点, 相切一个交点, 相交两个交点), 计算两圆相交(返回交点和个数), 过某点圆的切线(一条或两条), 两圆的切线(相离,内切,内含,外切)

#include <cstdio>

#include <iostream>

#include <cmath>

#include <cstdlib>

#include <cstring>

#include <vector>

#include <algorithm>

using namespace std;

#define PI  acos(-1)

const double eps = 1e-6;

struct Point

{

    double x, y;

    Point (double x = 0, double y = 0) : x(x), y(y) { } //构造函数, 方便代码书写

};

typedef Point myvector;

// 向量 + 向量 = 向量

myvector operator + (myvector A, myvector B) { return myvector(A.x + B.x, A.y + B.y); }

// 点 - 点 = 向量

myvector operator - (Point A, Point B) { return myvector(A.x - B.x, A.y - B.y); }

//向量 * 数 = 向量

myvector operator * (myvector A, double p) { return myvector(A.x * p, A.y * p); }

//向量/数 = 向量

myvector operator / (myvector A, double p) { return myvector(A.x / p, A.y / p); }

// 小于号

bool operator < (const Point & a, const Point & b)

{

    if (a.x == b.x) return a.y < b.y;

    return a.x < b.x;

}

//比较

int dcmp(double x)

{

    if(fabs (x) < eps) return 0;

    else return x < 0 ? -1 : 1;

}

// 恒等于号

bool operator ==  (const Point & a, const  Point  & b)

{

    return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0;

}

// 计算向量 A B 的点积, A*B = |A| * |B| * cosß

double Dot (myvector A, myvector B) { return A.x*B.x + A.y*B.y; }

// 计算向量 A 的长度

double Length (myvector A) { return sqrt (Dot(A, A)); }

// 计算向量 A,B 的夹角,是cos 有公式

double Angle (myvector A, myvector B)

{ return acos(Dot(A, B) / Length(A) / Length(B)); }

// 计算叉积,AxB = |A| * |B| * sinß, 得到的是与这两个向量垂直的向量

double Cross(myvector A, myvector B) { return A.x * B.y - A.y * B.x; }

double Area2(Point A, Point B, Point C) { return Cross (B - A, C - A); }

//计算两点距离

double DistancePoint(Point A, Point B) { return sqrt((A.x-B.x)*(A.x-B.x) + (A.y-B.y)*(A.y-B.y)); }

// 计算向量旋转后变成的另一个向量, rad 是弧度

//公式 x1 = x * cosß - y * sinß, y1 = x * sinß + y * cosß;

myvector Rotate(myvector A, double rad)

{

    return myvector(A.x * cos(rad) - A.y * sin(rad),

                  A.x * sin(rad) + A.y * cos(rad));

}

//计算向量的单位法线, 在调用前确保 A 不是零向量

myvector Normal(myvector A)

{

    double L = Length(A);

    return myvector(-A.y / L,  A.x / L);

}

//直线可以用直线上一点p1, 和方向向量V表示, 即 向量P = 点p1 + V;

//计算两直线的 交点 , 调用前确保两直线有交点

Point  GetLineInstersection(Point P, myvector v, Point Q, myvector w)

{

    myvector u = P - Q;

    double t = Cross(w, u) / Cross(v, w);

    return P + v * t;

}

//点到直线的距离

double DistanceToLine(Point P, Point A, Point B)

{

    myvector v1 = B - A, v2 = P - A;

    return fabs(Cross(v1, v2) / Length(v1));

}

// 点到线段的距离, 有两种可能, 一种点在线段上方, 这时候算垂直, 不在线段上方;

double DistanceToSegment(Point P, Point A, Point B)

{

    if( A == B) return Length(P-A); //如果线段是一个点

    myvector v1 = B - A, v2 = P - A, v3 = P - B;

    if(dcmp(Dot(v1, v2)) < 0)      return Length(v2);

    else if(dcmp(Dot(v1, v3)) > 0) return Length(v3);

    else return fabs(Cross(v1, v2)) / Length(v1);

}

//计算点在直线上投影的点

Point GetLineProjectoin(Point P, Point A, Point B)

{

    myvector v = B - A;

    return A + v * (Dot(v, P-A) / Dot(v, v));

}

struct Line

{

    Point v, p, e;

    Point point(double t)

    {

       return (p + v * t);

    }

};

struct Circle

{

    Point c;

    double r;

    Circle(Point _c=0,double _r=0):c(_c),r(_r){}

    Point point(double a)///根据圆心角算圆上的点

    {

        return Point(c.x+cos(a)*r,c.y+sin(a)*r);

    }

};

double angle(myvector V) {return atan2(V.y, V.x);}

//将直线AB沿法线方向平移d得到的直线EF,

myvector move_d(Point A, Point B, double d, Line& L)

{

    myvector C = B - A;

    C = C/Length(C);

    C = Rotate(C, PI/2);

    L.p = A + C * d;

    L.e = B + C * d;

    L.v = L.e - L.p;

    return (L.v);

}

//圆与直线的交点, 相离,没有交点, 相切一个交点, 相交两个交点

int getLineCircleInteresection(Line L, Circle C, double& t1, double& t2, vector<Point>& sol)

{

    //printf("##%f %f %f %f\n", L.p.x, L.p.y, L.e.x, L.e.y);

    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;

   // printf("delta = %.16f, esp = %.16f\n", delta, eps);

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

}

int getTangents(Point p, Circle C, myvector* v)

{

    myvector u = C.c - p;

    double dis = Length(u);

    if (dis < C.r)  return 0;

    else if (dcmp(dis - C.r) == 0)

    {

        v[0] = Rotate(u, PI / 2.0);

        return 1;

    }

    else

    {

        double ang = asin(C.r / dis);

        v[0] = Rotate(u, -ang);

        v[1] = Rotate(u, +ang);

        return 2;

    }

}

//两圆的切线条数, (1)重合,无数条,(2)两圆内含没有公共点没有切线,(3)两圆内切,有1条,

//(4)两圆相交有2条, (5)两圆外切,3条, (6)两圆相离,4条公切线

//返回切线条数, a[i],b[i]分别是第i条切线在圆A和B的切点

int getTangentsCircle(Circle A, Circle B, Point* a, Point* b)

{

    int cnt = 0;

    if(A.r < B.r)  //swap

    {

        Circle temp; Point *temp1 = NULL;

        A = temp; A = B; B = temp;

        a = temp1; a = b; b = temp1;

    }

    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++;

    }

    //有外共切线

    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(PI + base); cnt++;

    }

    else if(d2 > rsum * rsum)

    {

        double ang = acos(A.r + B.r) / 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;

}

bool cmp(Point A, Point B)

{

    if(A.x == B.x) return A.y < B.y;

    return A.x < B.x;

}

//三角形外接圆

void FUN1(double x1, double y1, double x2, double y2, double x3, double y3);

//三角形内切圆

void FUN2(double x1, double y1, double x2, double y2, double x3, double y3);

//过圆外某点切线的角度

void FUN3(double xc, double yc, double r, double xp, double yp);

//过某条直线外一点半径为r的圆

void FUN4(double xp, double yp, double x1, double y1, double x2, double y2, double r);

//和两条相交直线相切的半径为r的圆

void FUN5(double x1, double y1, double x2, double y2, double x3, double y3, double x4, double y4, double r);

//和两个相离的圆相切的圆

void FUN6(double x1, double y1, double r1, double x2, double y2, double r2, double r);

int main()

{

    //freopen("1.txt", "r", stdin);

    char s[50], s1[50] = "CircumscribedCircle", s2[50] = "InscribedCircle", s3[50] = "TangentLineThroughPoint",

    s4[50] = "CircleThroughAPointAndTangentToALineWithRadius", s5[50] = "CircleTangentToTwoLinesWithRadius",

    s6[50] = "CircleTangentToTwoDisjointCirclesWithRadius";

    while(~scanf("%s", s))

    {

       if(!strcmp(s,s1))

       {

           double x1, y1, x2, y2, x3, y3;

           cin >> x1 >> y1 >> x2 >> y2 >> x3 >> y3;

           FUN1(x1, y1, x2, y2, x3, y3);

       }

       else if(!strcmp(s, s2))

       {

           double x1, y1, x2, y2, x3, y3;

           cin >> x1 >> y1 >> x2 >> y2 >> x3 >> y3;

           FUN2(x1, y1, x2, y2, x3, y3);

       }

       else if(!strcmp(s,s3))

       {

           double xc, yc, r, xp, yp;

           cin >> xc >> yc >> r >> xp >> yp;

           FUN3(xc, yc, r, xp, yp);

       }

       else if(!strcmp(s,s4))

       {

           double xp, yp, x1, y1, x2, y2, r;

           cin >> xp >> yp >> x1 >> y1 >> x2 >> y2 >> r;

           FUN4( xp, yp, x1, y1, x2, y2, r);

       }

       else if(!strcmp(s,s5))

       {

           double x1, y1, x2, y2, x3, y3, x4, y4, r;

           cin >> x1 >> y1 >> x2 >> y2 >> x3 >> y3 >> x4 >> y4 >> r;

           FUN5(x1, y1, x2, y2, x3, y3, x4, y4, r);

       }

       else if(!strcmp(s,s6))

       {

           double x1, y1, r1, x2, y2, r2, r;

           cin >> x1 >> y1 >> r1 >> x2 >> y2 >> r2 >> r;

           FUN6(x1, y1, r1, x2, y2, r2, r);

       }

       getchar();

    }

    return 0;

}

void FUN1(double x1, double y1, double x2, double y2, double x3, double y3)

{//三角形外接圆

    Point A, B, C, D, E, F;

    myvector AB, BC, DE, DF;

    A.x = x1; A.y = y1; B.x = x2; B.y = y2; C.x = x3; C.y = y3;

E.x = (A.x + B.x)/2.0; E.y = (A.y + B.y)/2.0;

F.x = (B.x + C.x)/2.0; F.y = (B.y + C.y)/2.0;

    AB = B - A;       BC = C - B;

    DE = Normal(AB);  DF = Normal(BC);

    D = GetLineInstersection(E, DE, F, DF);

    double r = DistancePoint(B, D);

    printf("(%f,%f,%f)\n", D.x, D.y, r);

    return;

}

void FUN2(double x1, double y1, double x2, double y2, double x3, double y3)

{//三角形内切圆

    Point A, B, C;

    A.x = x1; A.y = y1; B.x = x2; B.y = y2; C.x = x3; C.y = y3;

    myvector v11 = B - A;

    myvector v12 = C - A;

    myvector v21 = A - B;

    myvector v22 = C - B;

    double ang1 = (angle(v11) + angle(v12)) / 2.0;

    double ang2 = (angle(v21) + angle(v22)) / 2.0;

    myvector vec1 = myvector(cos(ang1), sin(ang1));

    myvector vec2 = myvector(cos(ang2), sin(ang2));

    Point O = GetLineInstersection(A, vec1, B, vec2);

    double r = DistanceToLine(O, A, B);

    printf("(%f,%f,%f)\n", O.x, O.y, r);

}

void FUN3(double xc, double yc, double r, double xp, double yp)

{//过圆外某点切线的角度

    myvector vc[5];

    int len = getTangents(Point(xp, yp), Circle(Point(xc, yc), r), vc);

    double tmp[5];

    for (int i = 0; i < len; ++i)

    {

        double ang = angle(vc[i]);

        if (ang < 0) ang += PI;

        ang = fmod(ang, PI);

        tmp[i] = ang * 180 / PI;

    }

    sort(tmp, tmp + len);

    printf("[");

    for (int i = 0; i < len; ++i)

    {

        printf("%.6lf", tmp[i]);

        if (i != len - 1) printf(",");

    }

    printf("]\n");

    return;

}

void FUN4(double xp, double yp, double x1, double y1, double x2, double y2, double r)

{//过某条直线外一点半径为r的圆

    Line L1, L2;

    Point X, Y, P, Q, pp[10];

    double t1, t2;

    int k = 0;

    vector<Point>sol, sol2;

    X.x = x1; X.y = y1; Y.x = x2; Y.y = y2; P.x = xp; P.y = yp;

    Circle C(P, r);

    move_d(X, Y, -r, L1);

    move_d(X, Y, r, L2);

    int f = getLineCircleInteresection(L1, C, t1, t2, sol),

       f1 = getLineCircleInteresection(L2, C, t1, t2, sol2);

    printf("[");

    if(f == 1)

    {

       pp[k++] = sol[0];

    //  printf("(%f,%f)", sol[0].x, sol[0].y);

    }

    if(f == 2)

    {

       pp[k++] = sol[0]; pp[k++] = sol[1];

    //  printf("(%f,%f),(%f,%f)", sol[0].x, sol[0].y, sol[1].x, sol[1].y);

    }

    if(f1 == 1)

    {

       pp[k++] = sol2[0];

    //  if(f != 0) printf(",");

    //  printf("(%f,%f)", sol2[0].x, sol2[0].y);

    }

    if(f1 == 2)

    {

       pp[k++] = sol2[0];

       pp[k++] = sol2[1];

       //if(f != 0) printf(",");

       //printf("(%f,%f),(%f,%f)", sol2[0].x, sol2[0].y, sol2[1].x, sol2[1].y);

    }

    sort(pp,pp+k);

    for(int i=0;i<k;i++)

    {

       printf("(%f,%f)", pp[i].x, pp[i].y);

       if(i != k-1) printf(",");

    }

    printf("]\n");

    return;

}

void FUN5(double x1, double y1, double x2, double y2, double x3, double y3, double x4, double y4, double r)

{//和两条相交直线相切的半径为r的圆

    Line L1, L2;

    Point A, B, C, D, P, pp[10];

    int k = 0;

    A.x = x1; A.y = y1; B.x = x2; B.y = y2; C.x = x3; C.y = y3; D.x = x4; D.y = y4;

    move_d(A, B, r, L1); move_d(C, D, r, L2); P = GetLineInstersection(L1.p, L1.v, L2.p, L2.v);

    pp[k++] = P;

//  printf("[(%f,%f),", P.x, P.y);

    move_d(A, B, r, L1); move_d(C, D, -r, L2); P = GetLineInstersection(L1.p, L1.v, L2.p, L2.v);

    pp[k++] = P;

//  printf("(%f,%f),", P.x, P.y);

    move_d(A, B, -r, L1); move_d(C, D, r, L2); P = GetLineInstersection(L1.p, L1.v, L2.p, L2.v);

    pp[k++] = P;

//  printf("(%f,%f),", P.x, P.y);

    move_d(A, B, -r, L1); move_d(C, D, -r, L2); P = GetLineInstersection(L1.p, L1.v, L2.p, L2.v);

    pp[k++] = P;

//  printf("(%f,%f)]\n", P.x, P.y);

    sort(pp, pp+k);

    printf("[");

    for(int i=0;i<k;i++)

    {

       printf("(%f,%f)", pp[i].x, pp[i].y);

       if(i != k-1) printf(",");

    }

    printf("]\n");

    return;

}

void FUN6(double x1, double y1, double r1, double x2, double y2, double r2, double r)

{//和两个相离的圆相切的圆

    Point a, b, c, pp[10];

    int k = 0;

    a.x = x1; a.y = y1;b.x = x2; b.y = y2;

    Circle A(a,r1), B(b,r2);

    Circle C(a,r1+r), D(b,r2+r);

    vector<Point>sol;

    int t = getCircleCircleIntersection(C, D, sol);

    if(t == 1) pp[k++] = sol[0]; //printf("[(%f,%f)]\n", sol[0].x, sol[0].y);

    else if(t == 2)

    {

pp[k++] = sol[0]; pp[k++] = sol[1]; //printf("[(%f,%f),(%f,%f)]\n", sol[0].x, sol[0].y,sol[1].x, sol[1].y);

    }

    sort(pp, pp+k);

    printf("[");

    for(int i=0;i<k;i++)

    {

       printf("(%f,%f)", pp[i].x, pp[i].y);

       if(i != k-1) printf(",");

    }

    printf("]\n");    return;

}

 //Power by LJH

  

时间: 2024-10-29 02:48:45

计算几何模版的相关文章

三维计算几何模版

网上找了一个三维计算几何模版,完善了一下,使它能使用了... #include <cstdio> #include <cstring> #include <algorithm> using namespace std; /***********基础*************/ const double EPS=0.000001; typedef struct Point_3D { double x, y, z; Point_3D(double xx = 0, doubl

九野的计算几何模版

#include<cstdio> #include<iostream> #include<algorithm> #include<string.h> #include<math.h> using namespace std; #define point Point const double eps = 1e-8; const double PI = acos(-1.0); double ABS(double x){return x>0? x

寒假总结。。。

虽然刷了4章usaco,但跟没有刷一样,唯一有印象的就是判大素数模版,反着做背包的dp,反着做数位dp(逆向思维),对lis统计的判重,还有连通分量有了一点大概的印象,还有更加坚信了构造就是找规律,然后提高了一下代码能力,没了. 这个寒假最失败的一点就是把白天和夜晚倒过来过,白天睡觉,晚上刷题,这样好累,做题的效率好低..今天一定要把时间调整过来. 还有三天,我的寒假还没有结束!三天搞6个专题(其实是5个半,splay再写几道经典题就结束了,然后LCT应该很好学,和splay有联系),我的寒假才

HDU 3124 Moonmist(最近园对模板)

题目链接:http://acm.hdu.edu.cn/showproblem.php?pid=3124 寻找最近的两个圆的距离! 代码如下: //计算几何模版... //最近圆对 //HDU 3124 #include <cstdio> #include <cstring> #include <cstdlib> #include <cmath> #include <iostream> #include <algorithm> #inc

UVA12304 2D Geometry 110 in 1! 计算几何

计算几何: 堆几何模版就可以了.... Description Problem E 2D Geometry 110 in 1! This is a collection of 110 (in binary) 2D geometry problems. CircumscribedCircle x1 y1 x2 y2 x3 y3 Find out the circumscribed circle of triangle (x1,y1)-(x2,y2)-(x3,y3). These three poi

模板-计算几何

计算几何算法相关模版, 可能有错误, 省选前持续更正中 重要的不是模版内容, 而是提供算法的实现思路.

总结-计算几何

计算几何 做的题目很少, 而且模版也不熟. 所以还有这几天的时间, 把几个经典的算法弄熟弄懂, 模版能打出来就行了吧. 1. 凸包 2. 旋转卡壳 3. 半平面交 题目: 1. [codevs 1249] 多边形的面积 求多边形面积, 要理解叉积的意义. 2. [codevs 1298] 凸包周长 [codevs 3201] 奶牛代理商 XI 凸包周长 3. [codevs 1302] 小矮人(2002年CEOI中欧信息学奥赛) 旋转卡壳求对踵点 4. [codevs 3273] 两圆的交 这个

C++ 类模板三(类模版中的static关键字)

//类模版中的static关键字 #include<iostream> using namespace std; /* 类模板本质上是c++编译器根据类型参数创建了不同的类, c++编译器在利用类模板生成类的时候会为每个类生成一个static变量 那么对于类中的static关键字就非常好理解了 static关键字修饰的变量是属于类的 同一个类的对象共享类的static静态变量 类模板中的static修饰的变量数据类型必须是确定的 不可以是类型参数 因为静态变量在类对象之前初始化 这时候还没有通

!HDU 4380 三角屋内有奇数个宝藏的三角形有多少个-计算几何-(向量叉乘&amp;线段与点的关系&amp;暴力枚举)

题意:小明要买三座房子,这三个房子构成一个三角形,已知n个房子的坐标,任何三个房子都不在一条直线上,又已知有m个宝藏的坐标,问房子构成的三角形内有奇数个宝藏的三角形有多少个.数据范围:n(3~100),m(1~1000) 分析: 简单的计算几何.记住这题的做法. 三角形内的点的个数=上面的线段下面的点的个数 -- 下面两条线段下面的点的个数(或者下面一条线段减上面两条线段,看具体位置情况,所以直接取绝对值就好) n个点有n(n-1)/2条线段,不超过1W,枚举每条线段,再枚举每个宝藏的坐标(10