C++几何库.docx
- 文档编号:5362416
- 上传时间:2022-12-15
- 格式:DOCX
- 页数:22
- 大小:21.55KB
C++几何库.docx
《C++几何库.docx》由会员分享,可在线阅读,更多相关《C++几何库.docx(22页珍藏版)》请在冰豆网上搜索。
C++几何库
#include
#include"Geometry.h"
#include"math.h"
#include"Overlap.h"
#definemax(x,y)((x)>(y)?
(x):
(y))
#definemin(x,y)((x)<(y)?
(x):
(y))
/*常量定义*/
constdoubleINF=1E200;
constdoubleEP=1E-10;
constintMAXV=300;
constdoublePI=3.14159265;
//浮点误差的处理
intdblcmp(doubled)
{
if(fabs(d) return0; return(d>0)? 1: -1; } //返回两点之间欧氏距离 doubledist(d_POINTp1,d_POINTp2) { return(sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y))); } //判断两个点是否重合 boolequal_point(d_POINTp1,d_POINTp2) { return((abs(p1.x-p2.x) } /*(sp-op)*(ep-op)的叉积 r=multiply(sp,ep,op),得到(sp-op)*(ep-op)的叉积 r>0: sp在矢量opep的顺时针方向; r=0: opspep三点共线; r<0: sp在矢量opep的逆时针方向*/ doublemultiply(d_POINTsp,d_POINTep,d_POINTop) { return((sp.x-op.x)*(ep.y-op.y)-(ep.x-op.x)*(sp.y-op.y)); } doubleamultiply(d_POINTsp,d_POINTep,d_POINTop) { returnfabs((sp.x-op.x)*(ep.y-op.y)-(ep.x-op.x)*(sp.y-op.y)); } /*矢量(p1-op)和(p2-op)的点积 r=dotmultiply(p1,p2,op),得到矢量(p1-op)和(p2-op)的点积如果两个矢量都非零矢量 r<0: 两矢量夹角为锐角; r=0: 两矢量夹角为直角; r>0: 两矢量夹角为钝角*/ doubledotmultiply(d_POINTp1,d_POINTp2,d_POINTp0) { return((p1.x-p0.x)*(p2.x-p0.x)+(p1.y-p0.y)*(p2.y-p0.y)); } /*判断点p是否在线段l上 条件: (p在线段l所在的直线上)&&(点p在以线段l为对角线的矩形内)*/ boolonline(LINESEGl,d_POINTp) { 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))); } //返回点p以点o为圆心逆时针旋转alpha(单位: 弧度)后所在的位置 d_POINTrotate(d_POINTo,doublealpha,d_POINTp) { d_POINTtp; 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; returntp; } /*返回顶角在o点,起始边为os,终止边为oe的夹角(单位: 弧度) 角度小于pi,返回正值 角度大于pi,返回负值 可以用于求线段之间的夹角*/ doubleangle(d_POINTo,d_POINTs,d_POINTe) { doublecosfi,fi,norm; doubledsx=s.x-o.x; doubledsy=s.y-o.y; doubledex=e.x-o.x; doubledey=e.y-o.y; cosfi=dsx*dex+dsy*dey; norm=(dsx*dsx+dey*dey)*(dex*dex+dey*dey); cosfi/=sqrt(norm); if(cosfi>=1.0)return0; if(cosfi<=-1.0)return-3.1415926; fi=acos(cosfi); if(dsx*dey-dsy*dex>0)returnfi;//说明矢量os在矢量oe的顺时针方向 return-fi; } /*判断点C在线段AB所在的直线l上垂足P的与线段AB的关系 本函数是根据下面的公式写的,P是点C到线段AB所在直线的垂足 ACdotAB r=---------------------- ||AB||^2 (Cx-Ax)(Bx-Ax)+(Cy-Ay)(By-Ay) =---------------------------------------------------- L^2 rhasthefollowingmeaning: r=0P=A r=1P=B r<0PisonthebackwardextensionofAB r>1PisontheforwardextensionofAB 0 */ doublerelation(d_POINTc,LINESEGl) { LINESEGtl; tl.s=l.s; tl.e=c; returndotmultiply(tl.e,l.e,l.s)/(dist(l.s,l.e)*dist(l.s,l.e)); } //求点C到线段AB所在直线的垂足P d_POINTperpendicular(d_POINTp,LINESEGl) { doubler=relation(p,l); d_POINTtp; tp.x=l.s.x+r*(l.e.x-l.s.x); tp.y=l.s.y+r*(l.e.y-l.s.y); returntp; } /*求点p到线段l的最短距离 返回线段上距该点最近的点np注意: np是线段l上到点p最近的点,不一定是垂足*/ doubleptolinesegdist(d_POINTp,LINESEGl,d_POINT&np) { doubler=relation(p,l); if(r<0) { np=l.s; returndist(p,l.s); } if(r>1) { np=l.e; returndist(p,l.e); } np=perpendicular(p,l); returndist(p,np); } //求点p到线段l所在直线的距离 //请注意本函数与上个函数的区别 doubleptoldist(d_POINTp,LINESEGl) { returnabs(multiply(p,l.e,l.s))/dist(l.s,l.e); } /*计算点到折线集的最近距离,并返回最近点. 注意: 调用的是ptolineseg()函数*/ doubleptopointset(intvcount,d_POINTpointset[],d_POINTp,d_POINT&q) { inti; doublecd=double(INF),td; LINESEGl; d_POINTtq,cq; for(i=0;i { l.s=pointset[i]; l.e=pointset[i+1]; td=ptolinesegdist(p,l,tq); if(td { cd=td; cq=tq; } } q=cq; returncd; } /*判断圆是否在多边形内*/ boolCircleInsidePolygon(intvcount,d_POINTcenter,doubleradius,d_POINTpolygon[]) { d_POINTq; doubled; q.x=0; q.y=0; d=ptopointset(vcount,polygon,center,q); if(d elsereturnfalse; } /*返回两个矢量l1和l2的夹角的余弦(-1~1) 注意: 如果想从余弦求夹角的话,注意反余弦函数的值域是从0到pi*/ doublecosine(LINESEGl1,LINESEGl2) { doubled=((l1.e.x-l1.s.x)*(l2.e.x-l2.s.x)+(l1.e.y-l1.s.y)*(l2.e.y-l2.s.y)); return(d/(dist(l1.e,l1.s)*dist(l2.e,l2.s))); } //返回线段l1与l2之间的夹角 //单位: 弧度范围(-pi,pi) doublelsangle(LINESEGl1,LINESEGl2) { d_POINTo,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; returnangle(o,s,e); } //判断线段u和v相交(包括相交在端点处) boolintersect(LINESEGu,LINESEGv) { 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)); } //判断线段u和v相交(不包括双方的端点) boolintersect_A(LINESEGu,LINESEGv) { return((intersect(u,v))&& (! online(u,v.s))&& (! online(u,v.e))&& (! online(v,u.e))&& (! online(v,u.s))); } //判断线段v所在直线与线段u相交 //方法: 判断线段u是否跨立线段v boolintersect_l(LINESEGu,LINESEGv) { returnmultiply(u.s,v.e,v.s)*multiply(v.e,u.e,v.s)>=0; } //已知两点坐标p1和p2,求过这两点的直线解析方程: a*x+b*y+c=0(a>=0) LINEmakeLine(d_POINTp1,d_POINTp2) { LINEtl; intsign=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); returntl; } //根据直线解析方程返回直线的斜率k,水平线返回0,竖直线返回1e200 doubleslope(LINEl) { if(abs(l.a)<1e-20)return0; if(abs(l.b)<1e-20)returnINF; return-(l.a/l.b); } //返回直线的倾斜角alpha(0-pi) //注意: atan()返回的是-PI/2~PI/2 doublealpha(LINEl) { if(abs(l.a) if(abs(l.b) doublek=slope(l); if(k>0) returnatan(k); else returnPI+atan(k); } //求点p关于直线l的对称点 d_POINTsymmetry(LINEl,d_POINTp) { d_POINTtp; 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); returntp; } //如果两条直线l1(a1*x+b1*y+c1=0),l2(a2*x+b2*y+c2=0)相交,返回true,且返回交点p boollineInterSect(LINEl1,LINEl2,d_POINT&p)//是L1,L2 { doubled=l1.a*l2.b-l2.a*l1.b; if(abs(d) returnfalse; p.x=(l2.c*l1.b-l1.c*l2.b)/d; p.y=(l2.a*l1.c-l1.a*l2.c)/d; returntrue; } //如果线段l1和l2相交,返回true且交点由(inter)返回,否则返回false boolintersection(LINESEGl1,LINESEGl2,d_POINT&inter) { LINEll1,ll2; ll1=makeLine(l1.s,l1.e); ll2=makeLine(l2.s,l2.e); if(lineInterSect(ll1,ll2,inter)) returnonline(l1,inter); else returnfalse; } //如果无特别说明,输入多边形顶点要求按逆时针排列 //返回多边形面积(signed); //输入顶点按逆时针排列时,返回正值;否则返回负值 doublearea_of_polygon(intvcount,d_POINTpolygon[]) { inti; doubles; if(vcount<3) return0; s=polygon[0].y*(polygon[vcount-1].x-polygon[1].x); for(i=1;i s+=polygon[i].y*(polygon[(i-1)].x-polygon[(i+1)%vcount].x); returns/2; } //判断顶点是否按逆时针排列 //如果输入顶点按逆时针排列,返回true boolisconterclock(intvcount,d_POINTpolygon[]) { returnarea_of_polygon(vcount,polygon)>0; } /*射线法判断点q与多边形polygon的位置关系 要求polygon为简单多边形,顶点时针排列 如果点在多边形内: 返回0 如果点在多边形边上: 返回1 如果点在多边形外: 返回2*/ intinsidepolygon(intvcount,d_POINTPolygon[],d_POINTq) { intc=0,i; LINESEGl1,l2; l1.s=q;l1.e=q;l1.e.x=double(INF); //n=vcount; for(i=0;i { l2.s=Polygon[i]; l2.e=Polygon[(i+1)%vcount]; //doubleee=Polygon[(i+2)%vcount].x; //doubless=Polygon[(i+3)%vcount].y; if(online(l2,q)) return1; if(intersect_A(l1,l2)) c++;//相交且不在端点 if(online(l1,l2.e)&&! online(l1,l2.s)&&l2.e.y>l2.e.y) c++;//l2的一个端点在l1上且该端点是两端点中纵坐标较大的那个 if(! online(l1,l2.e)&&online(l1,l2.s)&&l2.e.y c++;//忽略平行边 } if(c%2==1) return0; else return2; } //判断点q在凸多边形polygon内 //点q是凸多边形polygon内[包括边上]时,返回true //注意: 多边形polygon一定要是凸多边形 boolInsideConvexPolygon(intvcount,d_POINTpolygon[],d_POINTq) { d_POINTp; LINESEGl; inti; p.x=0;p.y=0; for(i=0;i 多边形顶点平均值 { p.x+=polygon[i].x; p.y+=polygon[i].y; } p.x/=vcount; p.y/=vcount; for(i=0;i { l.s=polygon[i]; l.e=polygon[(i+1)%vcount]; if(multiply(p,l.e,l.s)*multiply(q,l.e,l.s)<0) /*点p和点q在边l的两侧,说明点q肯定在多边形外*/ returnfalse; } returntrue; } /*寻找凸包的graham扫描法 PointSet为输入的点集; ch为输出的凸包上的点集,按照逆时针方向排列; n为PointSet中的点的数目 len为输出的凸包上的点的个数*/ voidGraham_scan(d_POINTPointSet[],d_POINTch[],intn,int&len) { inti,j,k=0,top=2; d_POINTtmp; //选取PointSet中y坐标最小的点PointSet[k],如果这样的点有多个,则取最左边的一个 for(i=1;i if(PointSet[i].y &&(PointSet[i].x k=i; tmp=PointSet[0]; PointSet[0]=PointSet[k]; PointSet[k]=tmp;//现在PointSet中y坐标最小的点在PointSet[0] for(i=1;i 的按照距离PointSet[0]从近到远进行排序*/ { k=i; for(j=i+1;j if(multiply(PointSet[j],PointSet[k],PointSet[0])>0||//极角更小 (multiply(PointSet[j],PointSet[k],PointSet[0])==0)&&/*极角相等,距离更短*/dist(PointSet[0],PointSet[j]) k=j; tmp=PointSet[i]; PointSet[i]=PointSet[k]; PointSet[k]=tmp; } ch[0]=PointSet[0]; ch[1]=PointSet[1]; ch[2]=PointSet[2]; for(i=3;i { while(multiply(PointSet[i],ch[top],ch[top-1])>=0)top--; ch[++top]=PointSet[i]; } len=top+1; } //求凸多边形的重心,要求输入多边形按逆时针排序 d_POINTgravitycenter(intvcount,d_POINTpolygon[]) { d_POINTtp; doublex,y,s,x0,y0,cs,k; x=0;y=0;s=0; for(inti=1;i { x0=(polygon[0].x+polygon[i].x+polygon[i+1].x)/3; y0=(polygon[0].y+polygon[i].y+polygon[i+1].y)/3;//求当前三角形的重心 cs=multiply(polygon[i],polygon[i+1],polygon[0])/2; //三角形面积可以直接利用该公式求解 if(abs(s)<1e-20) { x=x0;y=y0;s+=cs;continue; } k=cs/s;//求面积比例 x=(x+k*x0)/(1+k); y=(y+k*y0)/(1+k); s+=cs; } tp.x=x; tp.y=y; returntp; } /*所谓凸多边形的直径,即凸多边形任两个顶点的最大距离。 下面的算法 仅耗时O(n),是一个优秀的算法。 输入必须是一个凸多边形,且顶点 必须按顺序(顺时针、逆时针均可)依次输入。 若输入不是凸多边形 而是一般点集,则要先求其凸包。 就是先求出所有跖对,然后求出每 个跖对的距离,取最大者。 点数要多于5个*/ voidDiameter(d_POINTch[],intn,double&dia) { intznum=0,i,j,k=1; intzd[MAXV][2]; doubletmp
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- C+ 几何