Graham-Scan法求凸包（代码copy自浙大模板）
#define eps 1e-8
#define zero(x) (((x)>0?(x):-(x))<eps)
struct Point{double x,y;};
Point pt[1001], stk[1001];
//计算cross product (P1-P0)x(P2-P0)
double Xmult(Point p1,Point p2,Point p0){
return (p1.x-p0.x)*(p2.y-p0.y)-(p2.x-p0.x)*(p1.y-p0.y);
}
//Graham_Scan算法顺时针构造包含所有共线点的凸包,O(nlogn)
Point p1,p2;
int graham_cp(const void* a,const void* b){
double ret=Xmult(*((Point*)a),*((Point*)b),p1);
return zero(ret)?(Xmult(*((Point*)a),*((Point*)b),p2)>0?1:-1):(ret>0?1:-1);
}
void _graham(int n,Point* pt,int& s,Point* ch){
int i,k=0;
for (p1=p2=pt[0],i=1;i<n;p2.x+=pt[i].x,p2.y+=pt[i].y,i++)
if (p1.y-pt[i].y>eps||(zero(p1.y-pt[i].y)&&p1.x>pt[i].x))
p1=pt[k=i];
p2.x/=n,p2.y/=n;
pt[k]=pt[0],pt[0]=p1;
qsort(pt+1,n-1,sizeof(Point),graham_cp);
for (ch[0]=pt[0],ch[1]=pt[1],ch[2]=pt[2],s=i=3;i<n;ch[s++]=pt[i++])
for (;s>2&&Xmult(ch[s-2],pt[i],ch[s-1])<-eps;s--);
}
//原始点集pt，凸包点集convex，返回凸包点集大小
int Graham_Scan(int n,Point* pt,Point* convex,int iscollinear=1,int isclockwise=1){
Point* temp=new Point[n];
int s,i;
_graham(n,pt,s,temp);
for (convex[0]=temp[0],n=1,i=(isclockwise?1:(s-1));isclockwise?(i<s):i;i+=(isclockwise?1:-1))
if (iscollinear||!zero(Xmult(temp[i-1],temp[i],temp[(i+1)%s])))
convex[n++]=temp[i];
delete []temp;
return n;
}
判断线段是否相交（包括规范相交），如果相交，求出交点，代码参考黑书。
#define Type int
#define fab(x) abs(x)
#define eps 1e-8
struct Point{
Type x, y;
};
int cmp( Type d ){
if( fab(d)< eps ) return 0;
return (d>0?1:-1);
}
int xyCmp( Type p, Type mini, Type maxi ){
return cmp(p-mini)*cmp(p-maxi);
}
int betweenCmp( Point a, Point b, Point c ){
if( fab(b.x-c.x)> fab(b.y-c.y) )
return xyCmp(a.x,min(b.x,c.x),max(b.x,c.x));
else
return xyCmp(a.y,min(b.y,c.y),max(b.y,c.y));
}
Type det( Type x1, Type y1, Type x2, Type y2 ){
return x1*y2-x2*y1;
}
Type cross( Point a, Point b, Point c ){
return det(b.x-a.x,b.y-a.y,c.x-a.x,c.y-a.y);
}
int SegCross( Point a, Point b, Point c, Point d, Point& p ){
Type s1, s2, s3, s4;
int d1= cmp(s1=cross(a,b,c));
int d2= cmp(s2=cross(a,b,d));
int d3= cmp(s3=cross(c,d,a));
int d4= cmp(s4=cross(c,d,b));
//规范相交求交点
if( (d1^d2)==-2 && (d3^d4)==-2 ){
p.x= ( c.x*s2-d.x*s1 )/ ( s2-s1 );
p.y= ( c.y*s2-d.y*s1 )/ ( s2-s1 );
return 1;
}
//判定非规范相交
if( d1==0 && betweenCmp(c,a,b)<=0 ||
d2==0 && betweenCmp(d,a,b)<=0 ||
d3==0 && betweenCmp(a,c,d)<=0 ||
d4==0 && betweenCmp(b,c,d)<=0 )
return 2;
return 0;
}
posted on 2009-06-26 13:26
古月残辉 阅读(683)
评论(1) 编辑 收藏 引用 所属分类:
计算几何