1、計算幾何幾何函數庫導引1.2.3.常量定義和包含文件基本數據結構精度控制 點的基本運算 1.2.3.4.5.6.7.平面上兩點之間距離判斷兩點是否重合 矢量叉乘矢量點乘判斷點是否段上求一點饒某點旋轉后的坐標求矢量夾角 線段及直線的基本運算 1.2.3.4.5.6.7.8.9.點與線段的關系求點到線段所在直線垂線的垂足點到線段的最近點點到線段所在直線的距離點到折線集的最近距離 判斷圓是否在多邊形內 求矢量夾角余弦求線段之間的夾角判斷線段是否相交判斷線段是否相交但不交在端點處求點關于某直線的對稱點判斷兩條直線是否相交及求直線交點判斷線段是否相交，如果相交返回交點 多邊形常用算法模塊 1.2.3.4

2、.5.7.8.9.判斷多邊形是否簡單多邊形檢查多邊形頂點的凸凹性 判斷多邊形是否凸多邊形 求多邊形面積判斷多邊形頂點的排列方向射判斷點是否在多邊形內判斷點是否在凸多邊形內尋找點集的 graham 算法10.尋找點集凸包的卷法凸包 MelkMan 算法的實現凸多邊形的直徑求凸多邊形的重心=導引/* 需要包含的頭文件 */#include /* 常量定義 */const double INF = 1E200; const double EP = 1E-10;constMAXV = 300;const doubl= 3.14159265;/* 基本幾何結構 */struct POdouble x;

3、double y;PO(double a=0, doub=0) x=a; y=b;struct LIPO PO LI LI;/ 直線的struct LINEGs;e;G(PO G() a, POb) s=a; e=b;方程 a*x+b*y+c=0為表示，約定 a= 0double a;doub;double c;LINE(double d1=1, double d2=-1, double d3=0) a=d1; b=d2; c=d3;/線段樹struct LINETREE/浮點誤差的處理dblcmp(double d)ibs(d)0) ?1 :-1 ;點的基本運算/ 返回兩點之間歐氏距離dou

4、ble dist(POp1,POp2)return( sqrt( (p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y) ) );/ 判斷兩個點是否重合bool equal_po(POp1,POp2)return ( (abs(p1.x-p2.x)EP)&(abs(p1.y-p2.y)0:sp 在矢量 op ep 的順時針方向；r=0：op sp ep 三點共線；r0: sp 在矢量 op ep 的逆時針方向 */double multiply(POsp,POep,POop)return(sp.x-op.x)*(ep.y-op.y) - (ep.x-op

5、.x)*(sp.y-op.y);double amultiply(POsp,POep,POop)return fabs(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:兩矢量夾角為鈍角 */double dotmultiply(POp1,POp2,POp0)return (p1.x-p0.x)*(p2.x-p0.x) + (p1.y-p0.y)*(p2.y-p0.y);/* 判斷點p

6、是否段 l 上條件：(p段 l 所在的直線上)& (點 p 在以線段 l 為對角線的矩形內) */bool online(LIG l,POp)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) =1.0 ) return 0;if (cosfi 0) return fi;/ 說明矢量 os 在矢量 oe 的順時針方向return -fi;線段及直線的基本運算/* 判斷點C段 AB 所在的直線 l 上垂足 P 的與線段 AB 的關系本函數是根據下面的公式寫的，

7、P 是點 C 到線段 AB 所在直線的垂足AC dot ABr =|AB|2(Cx-Ax)(Bx-Ax) + (Cy-Ay)(By-Ay)=L2r has the following meaning:r=0 r=1 r1 0r1P = A P = BP is on the backward exten P is on the forward extenP iserior to ABof AB of AB*/double relation(POc,LIG l)LItl.s=l.s; tl.e=c;G tl;return dotmultiply(tl.e,l.e,l.s)/(dist(l.s,l.

8、e)*dist(l.s,l.e);/ 求點 C 到線段 AB 所在直線的垂足 PPOpendicular(POp,LIG l)double r=relation(p,l);POtp;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;/* 求點 p 到線段 l 的最短距離返回線段上距該點最近的點 np注意：np 是線段 l 上到點 p 最近的點，不一定是垂足*/double ptoligdist(POp,LIG l,PO&np)double r=relation(p,l); if(r1)np=l.e;return d

9、ist(p,l.e);np=pendicular(p,l);return dist(p,np);/ 求點 p 到線段 l 所在直線的距離/請注意本函數與上個函數的區別double ptoldist(POp,LIG l)return abs(multiply(p,l.e,l.s)/dist(l.s,l.e);/* 計算點到折線集的最近距離,并返回最近點.注意：調用的是 ptolig()函數 */vcount, POdouble ptopoi;set(poset, POp, PO&q)double cd=double(INF),td;LIG l;POtq,cq;for(i=0;ivcount-1;

10、i+)l.s=po l.e=po td=ptoli if(tdcd)cd=td; cq=tq;seti;seti+1;gdist(p,l,tq);q=cq; return cd;/* 判斷圓是否在多邊形內*/bool CircleInsidePolygon(vcount,POcenter,double radius,POpolygon)POq;double d; q.x=0;q.y=0; d=ptoposet(vcount,polygon,center,q);if(dradius|fabs(d-radius)=min(v.s.x,v.e.x)&(max(v.s.x,v.e.x)=min(u.s

11、.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 相交（不包括雙方的端點）boolersect_A(LIG u,LIG v)return (ersect(u,v) &(!online(u,v.s)(!online(u,v.e)(!online(v,u.

12、e)(!online(v,u.s);& & &/ 判斷線段 v 所在直線與線段 u 相交方法：判斷線段 u 是否跨立線段 vboolersect_l(LIG u,LIG v)return multiply(u.s,v.e,v.s)*multiply(v.e,u.e,v.s)=0;/ 根據已知兩點坐標，求過這兩點的直線方程： a*x+b*y+c = 0(a = 0)LINE makeline(POLINE tl;sign = 1; tl.a=p2.y-p1.y; if(tl.a0)sign = -1; tl.a=sign*tl.a;p1,POp2)tl.b=sign*(p1.x-p2.x); t

13、l.c=sign*(p1.y*p2.x-p1.x*p2.y); return tl;/ 根據直線方程返回直線的斜率 k,水平線返回0,豎直線返回 1e200double slope(LINE l)if(abs(l.a) 1e-20)return 0; if(abs(l.b) 1e-20)return INF; return -(l.a/l.b);/ 返回直線的傾斜角 alpha ( 0 - pi)/ 注意：atan()返回的是 -PI/2 PI/2double alpha(LINE l)if(abs(l.a) EP)return 0; if(abs(l.b)0)return a);elsere

14、turn PI+a);/ 求點 p 關于直線 l 的對稱點POsymmetry(LINE l,POp)POtp;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;/ 如果兩條直線 l1(a1*x+b1*y+c1 = 0), l2(a2*x+b2*y+c2 = 0)相交，返回 true，且返回交點 pbool lineersect(LINE l1,

15、LINE l2,PO&p) / 是 L1，L2double d=l1.a*l2.b-l2.a*l1.b; if(abs(d)EP) / 不相交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;/ 如果線段 l1 和 l2 相交，返回 true 且交點由(er)返回，否則返回 falseboolersection(LIG l1,LIG l2,PO&er)LINE ll1,ll2; ll1=makeline(l1.s,l1.e); ll2=makeline(l2.s,l2.e)

16、;if(lineersect(ll1,ll2,er) return online(l1,er);else return false; 多邊形常用算法模塊如果無特別說明，輸入多邊形頂點要求按逆時針排列/ 返回多邊形面積(signed)；/ 輸入頂點按逆時針排列時，返回正值；否則返回負值double area_of_polygon(i; double s;if (vcount3) return 0;vcount,POpolygon)s=polygon0.y*(polygonvcount-1.x-polygon1.x); for (i=1;i0;/*射判斷點 q 與多邊形 polygon 的位置關系

17、要求 polygon 為簡單多邊形，頂點時針排列如果點在多邊形內： 返回 0如果點在多邊形邊上：返回 1如果點在多邊形外： 返回 2 */insidepolygon(POq)c=0,i,n;G l1,l2;LIl1.s=q; l1.e=q;l1.e.x=double(INF); n=vcount;for (i=0;il2.e.y)c+;/l2 的一個端點在 l1 上且該端點是兩端點中縱坐標較大的那個if(!online(l1,l2.e)& online(l1,l2.s) & l2.e.yl2.e.y) c+;/忽略平行邊if(c%2 = 1)return 0; elsereturn 2;/判斷

18、點 q 在凸多邊形 polygon 內/ 點 q 是凸多邊形 polygon 內包括邊上時，返回 true/ 注意：多邊形 polygon 一定要是凸多邊形bool InsideConvexPolygon(vcount,POpolygon,POq)PO LIp;G l;i;p.x=0; p.y=0;for(i=0;ivcount;i+) / 尋找一個肯定在多邊形 polygon 內的點 p：多邊形頂點平均值p.x+=polygoni.x; p.y+=polygoni.y;p.x /= vcount;p.y /= vcount;for(i=0;ivcount;i+)l.s=polygoni; l

19、.e=polygon(i+1)%vcount; if(multiply(p,l.e,l.s)*multiply(q,l.e,l.s)0)/* 點 p 和點 q 在邊 l 的兩側，說明點 q 肯定在多邊形外return false；return true；*/*尋找凸包的 graham 掃描法PoSet 為輸入的點集；ch 為輸出的凸包上的點集，按照逆時針方向排列;n 為 PoSet 中的點的數目len 為輸出的凸包上的點的個數 */void Graham_scan(POi,j,k=0,top=2;PoSet,POch,n,&len)POtmp;/ 選取 PoSet 中 y 坐標最小的點 PoS

20、etk，如果這樣的點有多個，則取最左邊的一個for(i=1;in;i+)if ( PoSeti.yPoSetk.y | (PoSeti.y=PoSetk.y)& (PoSeti.xPoSetk.x) ) k=i;tmp=PoSet0;Po PoSet0=PoSetk;Setk=tmp; / 現在 PoSet 中 y 坐標最小的點在 PoSet0for (i=1;in-1;i+) /* 對頂點按照相對 PoSet0的極角從小到大進行排序，極角相同的按照距離 Pok=i;Set0從近到遠進行排序 */for (j=i+1;j0 |/ 極角更小(multiply(PoSetj,PoSetk,PoSe

21、t0)=0)&/* 極角相等，距離更短 */dist(PoSet0,Po k=j;Setj)dist(PoSet0,PoSetk) )tmp=PoSeti;Po PoSeti=Po Setk=tmp;Setk;ch0=Po ch1=Po ch2=PoSet0;Set1;Set2;for (i=3;i=0) top-;ch+top=PoSeti;len=top+1;/ 卷法求點集凸殼，參數說明同 graham 算法void ConvexClosure(POtop=0,i,index,PoSet,POch,n,&len);double curmax,curcos,curdis;PO LItmp;G l1,l2;bool useMAXV; tmp=Po Set0; index=0;/ 選取 y 最小點，如果多于一個，則選取最左點for(i=1;in;i+)if(PoSeti.ytmp.y|PoSeti.y = tmp.y&PoSeti.xtmp.x)index=i;usei=false;tmp=PoSetindex;=index; useindex=true;index=-1; chtop+=tmp; tmp.x-=100; l1.s=tmp; l1.e=ch0;l2.s=ch0;while(index!=)curmax=-100; curdis=0;/ 選