圆的相关计算(刘汝佳版)

来源:互联网 发布:arp攻击软件 编辑:程序博客网 时间:2024/04/29 18:40
#include <iostream>#include <cstdio>#include <cstring>#include <algorithm>#include <cmath>using namespace std;/*Point的模板*/struct Point{double x,y;Point(double x=0,double y=0):x(x),y(y){}};typedef Point Vector;//向量+向量=向量,点+向量=点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); }bool operator < (const Point& a, const Point& b) {return a.x < b.x || (a.x == b.x && a.y < b.y);}//比较const double eps = 1e-10;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;}//基本计算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); } //有向面积//A向量逆时针旋转α rad //x'=xcosα-ysinα;//y'=xsinα+ycosα;Vector Rotate(Vector A, double rad){return Vector(A.x*cos(rad)-A.y*sin(rad),A.x*sin(rad)+A.y*cos(rad));}//A的单位法线,也就是逆时针90°,长度变为1,注意A要非零向量Vector Normal(Vector A){ double L=Length(A);return Vector(-A.y/L,A.x/L);}/*Point的模板*//*Circle的模板*/const double PI=acos(-1);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);}};struct Line{Vector v; //参数方程L=p+v*tPoint p;  Line(Point p,Vector v):v(v), p(p) {}Point point(double t){     //根据参数方程L=v+(p-v)*t中t求出线上的点return p+v*t;}};//求直线L与圆C的交点,两个交点参数方程值t1,t2,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));t2 = (-f+sqrt(delta)) / (2*e); sol.push_back(L.point(t2));return 2;}double Angle(Vector v) { return atan2(v.y,v.x); }   //求v的极角//求圆与圆的交点,函数返回交点数,sol保存交点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);     //求极角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;}//过点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) {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);return 2;}}//求两圆的公切线,返回切线条数,-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); }  //保证圆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=Angle(B.c-A.c);// 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;}//有外共切线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;}/*Circle的模板*/

0 0