基于二值距离变换的图像细化的代码实现

二值图像的细化是指将二值图像的拓扑结构计算出来，基于拓扑结构的特征提取及其他处理，计算量将得到大大降低。图像细化的算法很多，这里贴出本人的个人的一些想法，即利用二值图像的距离变换图来得到二值图像的拓扑结构。

二值图像的距离变换是指计算二值图像中非零像素点到最近的零像素点的距离，距离的定义很多，这里直接调用opencv里的开源函数：

distanceTransform（）;

该函数一共有四个参数，具体参考cv如下：

其中距离类型有如下几种：

enum
{
    CV_DIST_USER    =-1,  /**< User defined distance */
    CV_DIST_L1      =1,   /**< distance = |x1-x2| + |y1-y2| */
    CV_DIST_L2      =2,   /**< the simple euclidean distance */
    CV_DIST_C       =3,   /**< distance = max(|x1-x2|,|y1-y2|) */
    CV_DIST_L12     =4,   /**< L1-L2 metric: distance = 2(sqrt(1+x*x/2) - 1)) */
    CV_DIST_FAIR    =5,   /**< distance = c^2(|x|/c-log(1+|x|/c)), c = 1.3998 */
    CV_DIST_WELSCH  =6,   /**< distance = c^2/2(1-exp(-(x/c)^2)), c = 2.9846 */
    CV_DIST_HUBER   =7    /**< distance = |x|<c ? x^2/2 : c(|x|-c/2), c=1.345 */
};

本文采用简单的欧式距离，即 CV_DIST_L2;

测试原图和利用cv得到的距离图如下：

为了显示需要，已经将距离限制在（0-255）内，从图中可以看出，越靠近结构的“中心”，亮度越大，也就是距离越大。

而这些看起来较亮的地方正是我们所要提取的拓扑，那么如何提取这些亮点呢？直观的想法就是通过阈值，但是这种想法肯定会失败，因为在腿部最亮的点可能是100，但是身体大部分地方（非中间）也可能是100。

本文的思路是这样的：亮点是相对周围像素的（局部性），那么可能会想，只要在窗口内的所有像素都小于中心点像素，那么中心点像素不就是拓扑点了吗？很遗憾，拓扑点的领域内也是有可能有其他拓扑点的，而这些拓扑点的亮度是不确定的。既然窗口内的像素有比中心点大的，也有小的，那么如何定义什么情况下中心点才是拓扑点呢？看下图1-5：

其中粉色表示中心像素点，橘黄色表示领域内小于中心点的像素，深红色表示领域内大于中心点的像素，其中数字表示处于队里方向的像素值均小于中心点的标号，图1在对立方向（成180）没有均小于中心点的像素对，图2到图5的像素对分别为1,2,3,4；图1中，领域内对立方向中没有均小于中心点的像素对，说明其并不在“结构的中心位置”，所谓的中心位置至少有一个对立方向是均小于它的。符合要求的像素对越多，中心像素的顶点效应就越明显，即图2至图5依次增强。但在挑选拓扑点时，这样的像素对并不是越多越好，如图5，中心点虽然是最强顶点，但好像被孤立了，在其周围出现下一个连续拓扑点的概率较低，而图2、3、4中，深红色的像素都可以是候选的下一个拓扑点。所以在这里设定一个阈值，即领域内满足要求的像素对（对立方向的像素，且都比中心点的值小），可以根据需求，设置成不同大小，会得到不同的效果，理论上，值越小（大于0），拓扑结构越细，值越大，拓扑越粗，甚至消失。。。。

下面贴出细化的实现代码：

#include"head.h"

void TongThin(u08 *src, u08 *thin, s32 width, s32 height, s32 thinStructThreshold)
{   
    // 16 direction
    s32 sinTable[16] = {2, 2, 2, 1, 0, -1, -2, -2, -2, -2, -2, -1, 0, 1, 2, 2};
    s32 cosTable[16] = {0, 1, 2, 2, 2,  2,  2,  1,  0, -1, -2, -2, -2,-2,-2,-1};
    s32 surroundingPos[16]   = {0};
    u08 surroundingValue[16] = {0};

    for (s32 i = 2; i < height - 2; i++)
    {
        for (s32 j = 2; j < width - 2; j++)// 边界保护
        {
            s32 idx = i * width + j;

            if (0 == src[idx])
            {
                continue;
            }
            // 保存领域值
            for (s32 k = 0; k < 16; k++)
            {
                s32 step            = sinTable[k] * width + cosTable[k];
                surroundingPos[k]   = idx + step;
                surroundingValue[k] = src[surroundingPos[k]];
            }
            // 计算符合要求的像素对
            s32 count  = 0;
            for (s32 k = 0; k < 8; k++)
            {
                if (   (surroundingValue[k] - src[idx] < 0)
                    && (surroundingValue[8 + k] - src[idx] < 0))
                {
                    count++;
                }
            }
            // 比较
            if (count >= thinStructThreshold)
            {
                thin[idx] = 255;
            }

        }
    }
}

main函数如下：

#include "head.h"
#include "tongthin.h"
#include  "test.h"

void main()
{
        Mat frame = imread("C:\\Users\\Administrator\\Desktop\\image\\inputvideo\\6.jpg");
        Mat bin;
        Mat distransform(bin.size(), CV_32FC1);

        Mat frameGray;

        cvtColor(frame, frameGray, CV_BGR2GRAY);
        threshold(frameGray, bin,100, 255,1);
        // 距离变换
        distanceTransform(bin, distransform,CV_DIST_L2, 3);
        //得到最大距离
        f32 maxDis = 0;
        for ( s32 h = 0; h < bin.rows; h++)
        {
            for (s32 w = 0; w < bin.cols; w++)
            {
                f32 dis = distransform.at<float>(h, w);
                
                if (dis > maxDis)
                {
                    maxDis = dis;
                }
            }
        }
        // 规定化
        Mat disShow;
        disShow = Mat::zeros(bin.size(), CV_8UC1);
        for ( s32 h = 0; h < bin.rows; h++)
        {
            for (s32 w = 0; w < bin.cols; w++)
            {
                disShow.at<uchar>(h, w) = 
                     static_cast<uchar> (255 * distransform.at<float>(h, w) / maxDis);
            }
        }
        //细化
        Mat edges;
        edges = Mat::zeros(disShow.size(), CV_8UC1);
        double t = (double)getTickCount();
        s32    testTime     =   100;
        for (s32 i = 0; i < testTime; i++)
        {
              TongThin(disShow.data, edges.data, disShow.cols, disShow.rows, 2);
        }
        t = 1000*((double)getTickCount() - t)/getTickFrequency();
        cout << "Time is "<< t / testTime << " milliseconds."<< endl;

        // display
        imwrite("C:\\Users\\Administrator\\Desktop\\image\\inputvideo\\11.jpg",disShow);
        namedWindow("Extracted frame", 0);
        imshow("Extracted frame",frame);
        namedWindow("二值", 0);
        imshow("二值",bin);
        namedWindow("距离图", 0);
        imshow("距离图",disShow);
        namedWindow("细化", 0);
        imshow("细化",edges);

        setMouseCallback("距离图", on_mouseShowGray, &disShow);
        setMouseCallback("Extracted frame", on_mouseShowRgb, &frame);

        while (1)
        {
            imshow("距离图", disShow);
            imshow("Extracted frame",frame);
            waitKey(40);
        }
}

16邻域，8个像素对，分别设置阈值

thinStructThreshold为1,2,3,4,5,结果如下：

运行时间为5ms左右：

此文章只是对自己思考的一个记录，顺便分享，更多的思考和更佳的效果留给你们。