随笔-20  评论-16  文章-36  trackbacks-0
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 古月残辉 阅读(700) 评论(1)  编辑 收藏 引用 所属分类: 计算几何

评论:
# 黑书啊[未登录] 2010-08-25 10:36 | 代码疯子
用了很久,发现自己的模板居然对POJ1410过不了,错了。
现在想买黑书,但是想想就要退了,不买了。  回复  更多评论
  

只有注册用户登录后才能发表评论。
网站导航: 博客园   IT新闻   BlogJava   知识库   博问   管理