二维凸包convex hull之C++及OpenCV实现

原文地址：http://blog.csdn.net/xizhibei/article/details/7430385

 

打算接下来好好研究下算法（很明显，算法才是王道啊），然后尽量用直观的方式输出，于是用OpenCV画图成了不二首选，各位看官接下来看到一堆“XXX之C++及OpenCV实现”之类的标题就别见怪了～

另外还有个打算，看到自己写的东西被别人拿去占为己有，不爽，开始贴版权了^_^。

本文出处：http://blog.csdn.net/xizhibei

============================================================

今天就是二维凸包，算法导论中文版584页说的就是凸包，现在，让我们来实现它。

话说，凸包在很多地方有着重要的作用，如手势识别，需要识别出手的轮廓的凸包，二维或者三维区域的边界等等。

而对于凸包算法，其中最有名的莫过于Graham扫描算法，它的复杂度为nlog(n)，过程很优美，相信你看过运行过程你就会同样觉得了。

简单来说，这个算法的过程就是这样：

1.计算求得输入点x坐标最小（如果x相等，则比较y是不是最小）的点，作为第一个点

2.其它的点按照极角按顺时针排列，如果有共线的，取最近的那一个

3.依次将0，1，2三个坐标压入栈

4.对于剩余的点，i依次从3到最后，看次栈顶，栈顶，以及当前的i对应的点，如果是非左转动，那么弹栈

5.将i对应的点压栈，转到3

好了，现在开始实现：

下面是一些头文件/定义以及简单函数实现

其中值得一提的便是叉积，二维向量A和B的叉积就是A.x * B.y - A.y *B.x ，结果为正就表明A经过顺时针到达B（当然，小于180），反之则逆时针，凭借这种算法，就可以按极角排序了。

[cpp] view
plaincopyprint?

#include <cv.h>  

#include <cxcore.h>  

#include <highgui.h>  

  

#define POINT_NUM 100  

#define WIDTH 800  

#define HEIGHT 800  

#define zero 1e-12  

  

#define DIS(a,b) sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y))  

#define SGN(x) (fabs(x)<zero?0:(x>0?1:-1))  

#define CROSS(a,b,c) ((b.x-a.x)*(c.y-a.y)-(b.y-a.y)*(c.x-a.x))//叉积，用来判断旋转方向  

#define CMP(a,b) (a.x < b.x || SGN(a.x - b.x)==0 && a.y < b.y)//坐标的比较  

#define RAND (rand() % 100000 / 100000.0)//产生0-1之间的浮点数  

  

CvPoint p[POINT_NUM];  

CvPoint* hull_p = new CvPoint[POINT_NUM];//用来储存凸包上的点  

int hull_size = 0;  

IplImage* img;  

//下面简单实现了栈操作  

inline void push(CvPoint* S,CvPoint pt)  

{  

    S[hull_size++] = pt;  

}  

inline CvPoint pop(CvPoint* S)  

{  

    return S[--hull_size];  

}  

  

inline void swap(int x,int y)  

{  

    CvPoint pt = p[x];  

    p[x] = p[y];  

    p[y] = pt;  

}  

  

inline bool compare(CvPoint a,CvPoint b,CvPoint c)  

{  

    int tmp = SGN(CROSS(a,b,c));  

    if(tmp != 0)  

        return tmp > 0;  

    else//如果两点共线的话，就需要比较远近了  

        return DIS(a,b) < DIS(a,b);  

}  

//快排，极角的排序  

void sort(int l,int r)  

{  

    CvPoint tmp = p[(l + r) / 2];  

    int i = l;  

    int j = r;  

    do  

    {  

        while(compare(p[0],p[i],tmp))i++;  

        while(compare(p[0],tmp,p[j]))j--;  

        if(i <= j)  

        {  

            swap(i,j);  

            i++;  

            j--;  

        }  

    }while(i <=j);  

    if(i < r)sort(i,r);  

    if(j > l)sort(l,j);  

}  

然后就是重点了：
这里也跟上篇文章一样，为了看清运行过程，加上了画线函数，还有100微秒的延时，这样就能看清了

[cpp] view
plaincopyprint?

void draw_hull()  

{  

    int min = -1;  

    for(int j = 0;j < POINT_NUM;j++)//找出x坐标最小的，作为起始点  

    {  

        if(min == -1 || CMP(p[j],p[min]))  

            min = j;  

    }  

    if(min != 0)  

        swap(0,min);  

  

    sort(1,POINT_NUM - 1);//其他点排序  

  

    push(hull_p,p[0]);  

    push(hull_p,p[1]);  

    push(hull_p,p[2]);  

    for(int i = 3;i < POINT_NUM;i++)  

    {  

        while(CROSS(hull_p[hull_size - 2],hull_p[hull_size - 1],p[i]) < 0)//非左转  

        {  

            pop(hull_p);  

            cvLine(img,hull_p[hull_size - 1],p[i],cvScalar(255,0,255));//为了看清运行过程而加的  

            cvShowImage("Image",img);  

            cvWaitKey(100);  

        }  

        cvLine(img,hull_p[hull_size - 1],p[i],cvScalar(255,0,255));  

        push(hull_p,p[i]);  

    }  

  

    cvPolyLine(img,&hull_p,&hull_size,1,1,cvScalar(0,0,255),2);//最终画出凸包  

}  

显示图片：

[cpp] view
plaincopyprint?

void show_outcome()  

{  

    cvSet(img,cvScalar(255,255,255));  

    CvScalar color = cvScalar(0,0,0);  

    for(int i = 0;i < POINT_NUM;i++)//画出每个点，十字  

    {  

        int x = p[i].x;  

        int y = p[i].y;  

        cvLine(img,cvPoint(x - 5,y),cvPoint(x + 5,y),color,2);  

        cvLine(img,cvPoint(x,y - 5),cvPoint(x,y + 5),color,2);  

    }  

    draw_hull();  

    cvShowImage("Image",img);  

    cvWaitKey(0);  

}  

然后是主函数，这次的随机点生成换了个方式，点随机分布在左边的椭圆以及右边的圆上，不再象上次那样随机产生在矩形区域内了：

[cpp] view
plaincopyprint?

int main()  

{  

    img = cvCreateImage(cvSize(WIDTH,HEIGHT),IPL_DEPTH_8U,3);  

      

    srand((unsigned)time(NULL));  

    double phase = RAND * CV_PI * 2.0;  

    for (int i = 0; i < POINT_NUM / 2; i++) {  

        double r =  RAND * WIDTH / 4.0;  

        double theta = RAND * 1.5 * CV_PI + phase;  

        p[i] = cvPoint( WIDTH /4 + r * cos(theta), HEIGHT / 2 + 2 * r * sin(theta) );//椭圆  

    }  

    phase = RAND * CV_PI * 2.0;  

    for (int i = 0; i < POINT_NUM / 2; i++) {  

        double r =  RAND * WIDTH / 4.0;  

        double theta = RAND * 1.5 * CV_PI + phase;  

        p[i + POINT_NUM / 2] = cvPoint(WIDTH / 4 * 3 +  r * cos(theta), HEIGHT / 2 +  r * sin(theta));//圆  

    }  

  

    show_outcome();  

    delete [] hull_p;  

    return 0;  

}  

接下来看看效果：如果想看动态过程的画就去运行程序看吧～



参考：

http://www.cnblogs.com/Booble/archive/2011/03/10/1980089.html