纯C语言求点集的凸包程序（含边界提取）

来源：互联网发布：zencart 仿牌数据编辑：程序博客网时间：2024/06/14 06:22
求取点集凸包的数学原理为最简单的，在网上能够找到。
/******************************************************************************* * 文件名称  : GeoEnvelope.h * 当前版本  : 1.3* 作    者  : I_am_No3* 设计日期  : 2017年3月8日 * 内容摘要  : 离散点集的凸包轮廓********************************************************************************/#include <stdio.h>#include <stdlib.h>#include <math.h>#include <time.h>#define MaxNode 1000int stack[MaxNode];int top=0;int last=0;typedef struct POINT{    int x;    int y;}POINT;POINT point[MaxNode];POINT pointout[MaxNode];void swap(POINT point[],int i,int j){    POINT tmp;    tmp=point[i];    point[i]=point[j];    point[j]=tmp;}double multi(POINT p1,POINT p2,POINT p0) //叉乘{    return ((p1.x-p0.x)*(p2.y-p0.y)-(p1.y-p0.y)*(p2.x-p0.x));}double distence(POINT p1,POINT p2) //p1,p2的距离{    //return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y));    return sqrt ((double)((p1.x-p2.x)^2+(p1.y-p2.y)^2));}int cmp(const void *a,const void *b) {    POINT *c=(POINT *)a;    POINT *d=(POINT *)b;    double k=multi(*c,*d,point[0]);    if(k<0) return 1;    else if(k==0&&distence(*c,point[0])>=distence(*d,point[0])) return 1; //极角相同，比距离    else return -1;}void grahamscan(int n){top=0;    int i,u;    u=0;    for(i = 1;i<= n-1;i++) //找到最左下的点p0        if((point[i].y < point[u].y)||(point[i].y==point[u].y&&point[i].x < point[u].x))            u=i;    swap(point,0,u);    qsort(point+1,n-1,sizeof(point[0]),cmp); //按极角排序    for(i = 0;i <= 2;i++) stack[i] = i; //p0,p1,p2入栈    top=2;    for(i = 3;i <= n-1;i++) //最终凸包的各顶点的编号依次存在stack[]中。    {        while(multi(point[i],point[stack[top]],point[stack[top-1]])>=0) //弹出非左转的点     {            if(top==0)break;            top--;    }        top++;        stack[top] = i;    }}POINT* getline( POINT*pointmp , int n)//求取完凸包之后提取轮廓的程序，可以注释掉{last=0;int i,u,d,dx,dy,j,p,x,y,s1,s2;x=0;y=0;int interchange = 0;double derta,b;int temp;POINT tmp;for(i=0;i<n;i++)point[i]=pointmp[i];grahamscan(n);stack[top+1]=stack[0];for(p=0;p<=top;p++){u=stack[p]; d=stack[p+1]; dx=point[d].x-point[u].x; dy=point[d].y-point[u].y;if(dx!=0)derta=(double)dy/dx;elsederta=0;if(point[d].x>point[u].x){if(point[d].y>=point[u].y)   {if(abs(dx)>=abs(dy))for(j=0;j<abs(dx);j++){tmp.x=j+point[u].x; tmp.y=(int)(derta*j+point[u].y);pointout[last]=tmp;last++;}        else                    for(j=0;j<abs(dy);j++){tmp.y=j+point[u].y; tmp.x=(int)(j/derta+point[u].x);pointout[last]=tmp;last++;}}if(point[d].y<point[u].y)   {if(abs(dx)>=abs(dy))for(j=0;j<abs(dx);j++){tmp.x=j+point[u].x; tmp.y=(int)(derta*j+point[u].y);pointout[last]=tmp;last++;}        else                    for(j=0;j<abs(dy);j++){tmp.y=-j+point[u].y; tmp.x=(int)(-j/derta+point[u].x);pointout[last]=tmp;last++;}}}if(point[d].x<point[u].x){if(point[d].y>=point[u].y)   {if(abs(dx)>=abs(dy))for(j=0;j<abs(dx);j++){tmp.x=-j+point[u].x; tmp.y=(int)(-derta*j+point[u].y);pointout[last]=tmp;last++;}        else                    for(j=0;j<abs(dy);j++){tmp.y=j+point[u].y; tmp.x=(int)(j/derta+point[u].x);pointout[last]=tmp;last++;}}if(point[d].y<point[u].y)   {if(abs(dx)>=abs(dy))for(j=0;j<abs(dx);j++){tmp.x=-j+point[u].x; tmp.y=(int)(-derta*j+point[u].y);pointout[last]=tmp;last++;}        else                    for(j=0;j<abs(dy);j++){tmp.y=-j+point[u].y; tmp.x=(int)(-j/derta+point[u].x);pointout[last]=tmp;last++;}}  }if(point[d].x==point[u].x){if(point[d].y<point[u].y)for(j=0;j<abs(dy);j++){tmp.x=point[u].x; tmp.y=-j+point[u].y;pointout[last]=tmp;last++;}else                for(j=0;j<abs(dy);j++){tmp.x=point[u].x; tmp.y=j+point[u].y;pointout[last]=tmp;last++;}}}return pointout;}int getnumber(void){ return (last-1);}int main(){int starttime=clock();POINT point1[MaxNode];int i,j;j=0;for(i=0;i<=19;i++){for(j=0;j<=19;j++){point1[i*20+j].x = j;point1[i*20+j].y = i;}}getline( point1 ,400 );i=getnumber();int endtime=clock();//printf("time:%lf\n",last);printf("time:%lf\n",(endtime-starttime)/1000.0);system("Pause");return 0;}
1 0