计算几何模板大选

来源：互联网发布：淘宝卖家会员关系管理编辑：程序博客网时间：2024/04/30 08:49

1.两点之间距离
2.判断两点是否重合
3.叉积//可判断点在线段或直线的哪一侧
4.点积
5.判断点p是否在线段l上
6.返回点p以点o为圆心逆时针旋转alpha(单位：弧度)后所在的位置
7.返回顶角在o点，起始边为os，终止边为oe的夹角
8.判断点与线段的关系
9.求点C到线段AB所在直线的垂足 P
10.求点p到线段l的最短距离,并返回线段上距该点最近的点np
11.求点p到线段l所在直线的距离
12.计算点到折线集的最近距离,并返回最近点
13.判断圆是否在多边形内
14.返回两个矢量 l1和l2的夹角的余弦
15.返回线段l1与l2之间的夹角
16.如果线段u和v相交(包括相交在端点处)时,返回true
17.(线段u和v相交)&&(交点不是端点)时返回true
18.线段v所在直线与线段u相交时返回true
19.根据已知两点坐标，求过这两点的直线解析方程
20.根据直线解析方程返回直线的斜率k, 水平线返回 0, 竖直线返回 1e200
21.返回直线的倾斜角
22.求点p关于直线l的对称点
23.两直线相交返回true并返回交点p,不相交则返回false
24.如果线段l1和l2相交，返回true且交点由(inter)返回，否则返回false
25.返回值:点p在圆内(包括边界)时,返回true
26.三点确定一个圆,不能构成圆返回false
27.两圆位置关系
28.空间点到平面距离
29.点在直线同侧返回true
30.两个圆(已判断为相交或相切)的交点rp1,rp2
31.两相交圆公共面积
32.圆和直线(ax+by+c=0,a>=0)关系

33.三角形内切圆

34.过圆外一点的直线与圆的两个切点

<span style="font-size:14px;">#include<iostream>
#include<cmath>
#include<cstdio>
#define INF 1E200
using namespace std;
const double PI = 3.14159265;
const double eps=0.0000000001;
struct Point //点
{
double x,y;
Point (double a=0,double b=0):x(a),y(b) {}
};
struct Line_segment //线段
{
Point s,e;
Line_segment() {}
Line_segment(Point a,Point b):s(a),e(b) {}
};
struct Line //直线
{
double A,B,C;
Line(double A=1,double B=-1,double C=0):A(A),B(B),C(C) {}
};
inline double Max(double a,double b)
{
return a>b?a:b;
}
inline double Min(double a,double b)
{
return a<b?a:b;
}
//计算几何开始!
double Dist(Point a,Point b) //1.两点之间距离
{
return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}
bool equal_Point(Point a,Point b) //2.判断两点是否重合
{
return fabs(a.x-b.x)<eps&&fabs(a.y-b.y)<eps;
}
double multiply(Point sp,Point ep,Point op) //3.叉积//要判断点在直线哪一侧时,sp,ep为线段或者直线上两点,op为判断的点
{
/******************************************************************************
r=multiply(sp,ep,op),得到(sp-op) 和 ( ep-op)的叉积
r>0; ep在矢量op sp的逆时针方向//点在直线右边
r=0;op sp ep 三点共线；
r<0;ep在矢量op sp的顺时针方向//点在直线左边
*******************************************************************************/
return((sp.x-op.x)*(ep.y-op.y)-(ep.x-op.x)*(sp.y-op.y));
}
double dotmultiply(Point p1,Point p2,Point p0) //4.点积
{
/******************************************************************************
r=dotmultiply(p1,p2,op),得到矢量(p1-op)和(p2-op)的点积,如果两个矢量都非零矢量
r<0:两矢量夹角为钝角;
r=0:两矢量夹角为直角;
r>0:两矢量夹角为锐角;
*******************************************************************************/
return (p1.x-p0.x)*(p2.x-p0.x)+(p1.y-p0.y)*(p2.y-p0.y);
}
bool online(Line_segment l,Point p)//5.判断点p是否在线段l上 +3
{
/******************************************************************************
判断点p是否在线段l上
条件：(p在线段l所在的直线上) && (点p在以线段l为对角线的矩形内)
*******************************************************************************/
return( (multiply(l.e,p,l.s)==0) &&( ( (p.x-l.s.x)*(p.x-l.e.x)<=0 )&&( (p.y-l.s.y)*(p.y-l.e.y)<=0 ) ) );
}
Point Rotate(Point o,double alpha,Point p) // 6.返回点p以点o为圆心逆时针旋转alpha(单位：弧度)后所在的位置
{
Point tp;
p.x-=o.x;
p.y-=o.y;
tp.x=p.x*cos(alpha)-p.y*sin(alpha)+o.x;
tp.y=p.y*cos(alpha)+p.x*sin(alpha)+o.y;
return tp;
}
double angle(Point o,Point s,Point e) //7.返回顶角在o点，起始边为os，终止边为oe的夹角(单位：弧度)
{
/******************************************************************************
返回顶角在o点，起始边为os，终止边为oe的夹角(单位：弧度)
角度小于pi，返回正值
角度大于pi，返回负值
可以用于求线段之间的夹角
原理：
r = dotmultiply(s,e,o) / (dist(o,s)*dist(o,e))
r'= multiply(s,e,o)
r >= 1 angle = 0;
r <= -1 angle = -PI
-1<r<1 && r'>0 angle = arccos(r)
-1<r<1 && r'<=0 angle = -arccos(r)
********************************************************************************/
double cosfi,fi,norm;
double dsx = s.x - o.x;
double dsy = s.y - o.y;
double dex = e.x - o.x;
double dey = e.y - o.y;
cosfi=dsx*dex+dsy*dey;
norm=(dsx*dsx+dsy*dsy)*(dex*dex+dey*dey);
cosfi /= sqrt( norm );
if (cosfi >= 1.0 ) return 0;
if (cosfi <= -1.0 ) return -3.1415926;
fi=acos(cosfi);
if (dsx*dey-dsy*dex>0) return fi; // 说明矢量os 在矢量 oe的顺时针方向
return -fi;
}
double relation(Point p,Line_segment l) //8.判断点与线段的关系
{
/******************************************************************************
判断点与线段的关系,用途很广泛
本函数是根据下面的公式写的，P是点C到线段AB所在直线的垂足
AC dot AB
r = ---------
||AB||^2
(Cx-Ax)(Bx-Ax) + (Cy-Ay)(By-Ay)
= -------------------------------
L^2
r has the following meaning:
r=0 P = A
r=1 P = B
r<0 P is on the backward extension of AB
r>1 P is on the forward extension of AB
0<r<1 P is interior to AB
********************************************************************************/
Line_segment tl;
tl.s=l.s;
tl.e=p;
return dotmultiply(tl.e,l.e,l.s)/(Dist(l.s,l.e)*Dist(l.s,l.e));
}
Point perpendicular(Point p,Line_segment l) //9.求点C到线段AB所在直线的垂足 P
{
/******************************************************************************
求点C到线段AB所在直线的垂足 P
*******************************************************************************/
double r=relation(p,l);
Point tp;
tp.x=l.s.x+r*(l.e.x-l.s.x);
tp.y=l.s.y+r*(l.e.y-l.s.y);
return tp;
}
double ptolinesegdist(Point p,Line_segment l,Point &np) //10.求点p到线段l的最短距离,并返回线段上距该点最近的点np
{
/******************************************************************************
求点p到线段l的最短距离,并返回线段上距该点最近的点np
注意：np是线段l上到点p最近的点，不一定是垂足
*******************************************************************************/
double r=relation(p,l);
if(r<0)
{
np=l.s;
return Dist(p,l.s);
}
if(r>1)
{
np=l.e;
return Dist(p,l.e);
}
np=perpendicular(p,l);
return Dist(p,np);
}
double ptoldist(Point p,Line_segment l) // 11.求点p到线段l所在直线的距离,请注意本函数与上个函数的区别
{
return abs(multiply(p,l.e,l.s))/Dist(l.s,l.e);
}
double ptopointset(int vcount,Point pointset[],Point p,Point &q) //12.计算点到折线集的最近距离,并返回最近点.
{
/******************************************************************************
计算点到折线集的最近距离,并返回最近点.
注意：调用的是ptolineseg()函数
*******************************************************************************/
int i;
double cd=double(INF),td;
Line_segment l;
Point tq,cq;
for(i=0; i<vcount-1; i++)
{
l.s=pointset[i];
l.e=pointset[i+1];
td=ptolinesegdist(p,l,tq);
if(td<cd)
{
cd=td;
cq=tq;
}
}
q=cq;
return cd;
}
bool CircleInsidePolygon(int vcount,Point center,double radius,Point polygon[]) //13.判断圆是否在多边形内.ptolineseg()函数的应用
{
Point q;
double d;
q.x=0;
q.y=0;
d=ptopointset(vcount,polygon,center,q);
if(d<radius||fabs(d-radius)<eps)
return true;
else
return false;
}
double cosine(Line_segment l1,Line_segment l2) //14.返回两个矢量 l1和l2的夹角的余弦
{
/******************************************************************************
返回两个矢量l1和l2的夹角的余弦(-1 --- 1)注意：如果想从余弦求夹角的话，注意反余弦函数的定义域是从 0到pi
*******************************************************************************/
return (((l1.e.x-l1.s.x)*(l2.e.x-l2.s.x) +(l1.e.y-l1.s.y)*(l2.e.y-l2.s.y))/(Dist(l1.e,l1.s)*Dist(l2.e,l2.s)));
}
double lsangle(Line_segment l1,Line_segment l2) // 15.返回线段l1与l2之间的夹角单位:弧度范围(-pi，pi)
{
Point o,s,e;
o.x=o.y=0;
s.x=l1.e.x-l1.s.x;
s.y=l1.e.y-l1.s.y;
e.x=l2.e.x-l2.s.x;
e.y=l2.e.y-l2.s.y;
return angle(o,s,e);
}
bool intersect(Line_segment u,Line_segment v) // 16.如果线段u和v相交(包括相交在端点处)时,返回true
{
/******************************************************************************
如果线段u和v相交(包括相交在端点处)时,返回true
判断P1P2跨立Q1Q2的依据是: ( P1 - Q1 ) x ( Q2 - Q1 ) x ( Q2 - Q1 ) x ( P2 - Q1 ) >= 0
判断Q1Q2跨立P1P2的依据是: ( Q1 - P1 ) x ( P2 - P1 ) x ( P2 - P1 ) x ( Q2 - P1 ) >= 0
*******************************************************************************/
return((Max(u.s.x,u.e.x)>=Min(v.s.x,v.e.x))&& //排斥实验
(Max(v.s.x,v.e.x)>=Min(u.s.x,u.e.x))&&
(Max(u.s.y,u.e.y)>=Min(v.s.y,v.e.y))&&
(Max(v.s.y,v.e.y)>=Min(u.s.y,u.e.y))&&
(multiply(v.s,u.e,u.s)*multiply(u.e,v.e,u.s)>=0)&& //跨立实验
(multiply(u.s,v.e,v.s)*multiply(v.e,u.e,v.s)>=0));
}
bool intersect_A(Line_segment u,Line_segment v) // 17.(线段u和v相交)&&(交点不是端点)时返回true
{
return ((intersect(u,v))&&
(!online(u,v.s))&&
(!online(u,v.e))&&
(!online(v,u.e))&&
(!online(v,u.s)));
}
bool intersect_l(Line_segment u,Line_segment v)// 18.线段v所在直线与线段u相交时返回true
{
/******************************************************************************
线段v所在直线与线段u相交时返回true;方法;判断线段u是否跨立线段v
*******************************************************************************/
return multiply(u.s,v.e,v.s)*multiply(v.e,u.e,v.s)>=0;
}
Line makeline(Point p1,Point p2) // 19.根据已知两点坐标，求过这两点的直线解析方程
{
/******************************************************************************
根据已知两点坐标，求过这两点的直线解析方程Ax+By+C=0 (A>=0)
*******************************************************************************/
Line tl;
int sign = 1;
tl.A=p2.y-p1.y;
if(tl.A<0)
{
sign = -1;
tl.A=sign*tl.A;
}
tl.B=sign*(p1.x-p2.x);
tl.C=sign*(p1.y*p2.x-p1.x*p2.y);
return tl;
}
double slope(Line l) // 20.根据直线解析方程返回直线的斜率k, 水平线返回 0, 竖直线返回 1e200
{
if(abs(l.A) < 1e-20)
return 0;
if(abs(l.B) < 1e-20)
return INF;
return -(l.A/l.B);
}
double alpha(Line l) //21. 返回直线的倾斜角 alpha ( 0 - pi)
{
if(abs(l.A)< eps)
return 0;
if(abs(l.B)< eps)
return PI/2;
double k=slope(l);
if(k>0)
return atan(k);
else
return PI+atan(k);
}
Point symmetry(Line l,Point p) //22. 求点p关于直线l的对称点
{
Point tp;
tp.x=((l.B*l.B-l.A*l.A)*p.x-2*l.A*l.B*p.y-2*l.A*l.C)/(l.A*l.A+l.B*l.B);
tp.y=((l.A*l.A-l.B*l.B)*p.y-2*l.A*l.B*p.x-2*l.B*l.C)/(l.A*l.A+l.B*l.B);
return tp;
}
bool lineintersect(Line l1,Line l2,Point &p) // 23.两直线相交返回true并返回交点p,不相交则返回false
{
double d=l1.A*l2.B-l2.A*l1.B;
if(abs(d)<eps) // 不相交
return false;
p.x = (l2.C*l1.B-l1.C*l2.B)/d;
p.y = (l2.A*l1.C-l1.A*l2.C)/d;
return true;
}
bool intersection(Line_segment l1,Line_segment l2,Point &inter) // 24.如果线段l1和l2相交，返回true且交点由(inter)返回，否则返回false
{
Line ll1,ll2;
ll1=makeline(l1.s,l1.e);
ll2=makeline(l2.s,l2.e);
if(lineintersect(ll1,ll2,inter))
return online(l1,inter) && online(l2,inter);
else
return false;
}
/*******************************************************************************/
bool point_in_circle(Point o,double r,Point p) //25. 返回值:点p在圆内(包括边界)时,返回true
{
/******************************************************************************
参数o为圆心,r为半径,p为判断的点
返回值:点p在圆内(包括边界)时,返回true
*******************************************************************************/
double d2=(p.x-o.x)*(p.x-o.x)+(p.y-o.y)*(p.y-o.y);
double r2=r*r;
return d2<r2||abs(d2-r2)<eps;
}
bool cocircle(Point p1,Point p2,Point p3,Point &q,double &r) //26.三点确定一个圆,不能构成圆返回false
{
/******************************************************************************
用途 :求不共线的三点确定一个圆
输入 :三个点p1,p2,p3
返回值 :如果三点共线，返回false；反之，返回true。圆心由q返回，半径由r返回
*******************************************************************************/
double x12=p2.x-p1.x;
double y12=p2.y-p1.y;
double x13=p3.x-p1.x;
double y13=p3.y-p1.y;
double z2=x12*(p1.x+p2.x)+y12*(p1.y+p2.y);
double z3=x13*(p1.x+p3.x)+y13*(p1.y+p3.y);
double d=2.0*(x12*(p3.y-p2.y)-y12*(p3.x-p2.x));
if(abs(d)<eps) //共线，圆不存在
return false;
q.x=(y13*z2-y12*z3)/d;
q.y=(x12*z3-x13*z2)/d;
r=Dist(p1,q);
return true;
}
int CircleRelation(Point p1, double r1, Point p2, double r2) //27.两圆位置关系
{
/******************************************************************************
相离:return 1
外切:return 2
相交:return 3
内切:return 4
内含:return 5
*******************************************************************************/
double d = sqrt( (p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y) );
if( fabs(d-r1-r2) < eps ) // 必须保证前两个if先被判定！
return 2;
if( fabs(d-fabs(r1-r2)) < eps )
return 4;
if( d > r1+r2 )
return 1;
if( d < fabs(r1-r2) )
return 5;
if( fabs(r1-r2) < d && d < r1+r2 )
return 3;
return 0; // indicate an error!未知错误
}
double P2planeDist(double x, double y, double z, double a, double b, double c, double d) //28.空间点到平面距离
{
/******************************************************************************
空间点到平面的距离,平面用一般式表示ax+by+cz+d=0
*******************************************************************************/
return fabs(a*x+b*y+c*z+d) / sqrt(a*a+b*b+c*c);
}
bool SameSide(Point p1, Point p2, Line line) //29.点在直线同侧返回true
{
return (line.A * p1.x + line.B * p1.y + line.C) *
(line.A * p2.x + line.B * p2.y + line.C) > 0;
}
void c2point(Point p1,double r1,Point p2,double r2,Point &rp1,Point &rp2) //30.两个圆(已判断为相交或相切)的交点rp1,rp2
{
double a,b,r;
a=p2.x-p1.x;
b=p2.y-p1.y;
r=(a*a+b*b+r1*r1-r2*r2)/2;
if(a==0&&b!=0)
{
rp1.y=rp2.y=r/b;
rp1.x=sqrt(r1*r1-rp1.y*rp1.y);
rp2.x=-rp1.x;
}
else if(a!=0&&b==0)
{
rp1.x=rp2.x=r/a;
rp1.y=sqrt(r1*r1-rp1.x*rp2.x);
rp2.y=-rp1.y;
}
else if(a!=0&&b!=0)
{
double delta;
delta=b*b*r*r-(a*a+b*b)*(r*r-r1*r1*a*a);
rp1.y=(b*r+sqrt(delta))/(a*a+b*b);
rp2.y=(b*r-sqrt(delta))/(a*a+b*b);
rp1.x=(r-b*rp1.y)/a;
rp2.x=(r-b*rp2.y)/a;
}
rp1.x+=p1.x;
rp1.y+=p1.y;
rp2.x+=p1.x;
rp2.y+=p1.y;
}
double c2area(Point p1,double r1,Point p2,double r2) //31.两相交圆公共面积 +30
{
double TEMP;
Point rp1,rp2,rp;
c2point(p1,r1,p2,r2,rp1,rp2);
if(r1>r2) //保证r2>r1
{
rp=p1;
p1=p2;
p2=rp;
TEMP=r1;
r1=r2;
r2=TEMP;
}
double a,b,rr;
a=p1.x-p2.x;
b=p1.y-p2.y;
rr=sqrt(a*a+b*b);
double dx1,dy1,dx2,dy2;
double sita1,sita2;
dx1=rp1.x-p1.x;
dy1=rp1.y-p1.y;
dx2=rp2.x-p1.x;
dy2=rp2.y-p1.y;
sita1=acos((dx1*dx2+dy1*dy2)/r1/r1);
dx1=rp1.x-p2.x;
dy1=rp1.y-p2.y;
dx2=rp2.x-p2.x;
dy2=rp2.y-p2.y;
sita2=acos((dx1*dx2+dy1*dy2)/r2/r2);
double s=0;
if(rr<r2)//相交弧为优弧
s=r1*r1*(PI-sita1/2+sin(sita1)/2)+r2*r2*(sita2-sin(sita2))/2;
else//相交弧为劣弧
s=(r1*r1*(sita1-sin(sita1))+r2*r2*(sita2-sin(sita2)))/2;
return s;
}
\
int clpoint(Point p,double r,double a,double b,double c,Point &rp1,Point &rp2) //32.圆和直线(ax+by+c=0,a>=0)关系
{
/******************************************************************************
相离 return 0
相切 return 1
相交 return 2
*******************************************************************************/
int res=0;
c=c+a*p.x+b*p.y;
double tmp;
if(a==0&&b!=0)
{
tmp=-c/b;
if(r*r<tmp*tmp)
res=0;
else if(r*r==tmp*tmp)
{
res=1;
rp1.y=tmp;
rp1.x=0;
}
else
{
res=2;
rp1.y=rp2.y=tmp;
rp1.x=sqrt(r*r-tmp*tmp);
rp2.x=-rp1.x;
}
}
else if(a!=0&&b==0)
{
tmp=-c/a;
if(r*r<tmp*tmp)
res=0;
else if(r*r==tmp*tmp)
{
res=1;
rp1.x=tmp;
rp1.y=0;
}
else
{
res=2;
rp1.x=rp2.x=tmp;
rp1.y=sqrt(r*r-tmp*tmp);
rp2.y=-rp1.y;
}
}
else if(a!=0&&b!=0)
{
double delta;
delta=b*b*c*c-(a*a+b*b)*(c*c-a*a*r*r);
if(delta<0)
res=0;
else if(delta==0)
{
res=1;
rp1.y=-b*c/(a*a+b*b);
rp1.x=(-c-b*rp1.y)/a;
}
else
{
res=2;
rp1.y=(-b*c+sqrt(delta))/(a*a+b*b);
rp2.y=(-b*c-sqrt(delta))/(a*a+b*b);
rp1.x=(-c-b*rp1.y)/a;
rp2.x=(-c-b*rp2.y)/a;
}
}
rp1.x+=p.x;
rp1.y+=p.y;
rp2.x+=p.x;
rp2.y+=p.y;
return res;
}
void incircle(Point p1,Point p2,Point p3,Point &rp,double &r) // 33.三角形内切圆
{
double dx31,dy31,dx21,dy21,d31,d21,a1,b1,c1;
dx31=p3.x-p1.x;
dy31=p3.y-p1.y;
dx21=p2.x-p1.x;
dy21=p2.y-p1.y;
d31=sqrt(dx31*dx31+dy31*dy31);
d21=sqrt(dx21*dx21+dy21*dy21);
a1=dx31*d21-dx21*d31;
b1=dy31*d21-dy21*d31;
c1=a1*p1.x+b1*p1.y;
double dx32,dy32,dx12,dy12,d32,d12,a2,b2,c2;
dx32=p3.x-p2.x;
dy32=p3.y-p2.y;
dx12=-dx21;
dy12=-dy21;
d32=sqrt(dx32*dx32+dy32*dy32);
d12=d21;
a2=dx12*d32-dx32*d12;
b2=dy12*d32-dy32*d12;
c2=a2*p2.x+b2*p2.y;
rp.x=(c1*b2-c2*b1)/(a1*b2-a2*b1);
rp.y=(c2*a1-c1*a2)/(a1*b2-a2*b1);
r=fabs(dy21*rp.x-dx21*rp.y+dx21*p1.y-dy21*p1.x)/d21;
}
void cutpoint(Point p,double r,Point sp,Point &rp1,Point &rp2) //34.过圆外一点的直线与圆的两个切点(p为圆心,r为圆半径,点sp为圆外一点)
{
Point p2;
p2.x=(p.x+sp.x)/2;
p2.y=(p.y+sp.y)/2;
double dx2,dy2,r2;
dx2=p2.x-p.x;
dy2=p2.y-p.y;
r2=sqrt(dx2*dx2+dy2*dy2);
c2point(p,r,p2,r2,rp1,rp2);
}
int main()
{
double a,b,c,d;
Point m(1,1),n(2,2);
Line_segment w(m,n);
while(~scanf("%lf %lf",&a,&b))
{
Point A(a,b);
if(online(w,A)) printf("yes\n");
else printf("n0\n");
}
return 0;
}</span>

0 0