Apriltag源码学习

写这篇文章:记录一下项目中用到的AprilTag码部分内容。

本来以为是简单的应用一下,结果发现环境四边形多时,检测时间真的太久了,根本满足不了要求,只能硬着头皮看源码。这篇文章就是看源码后的一些总结。

在测试中发现基本上大部分时间都花在了算法的quad检测部分,所以我主要是学习了前半部分代码(其实是因为后面的操作不感兴趣🙂)

开始代码部分(看我这部分你需要打开Apirltag代码一起跟着走)

首先介绍一下AprilTag中我认为比较重要的结构体变量。

这个是AprilTag最基本的单元 类似于 vector<类型> var; var.reserve(alloc);
typedef struct zarray zarray_t;
struct zarray
{
    size_t el_sz;  // 记录每一个元素的大小

    int size; // 当前有多少的元素 可以用来计算剩余的内存空间 具体 alloc - size * int(el_sz)
    int alloc; // 我们分配了多少的内存大小 
    char *data; // 以字节为单位进行数据存储 (说白了就是容器, 没有啥具体含义) 1Byte 1Byte ... 
};

zarray_t的一些重要方法

static inline zarray_t *zarray_create(size_t el_sz) //  size_t el_sz 表示即将存储的元素大小
                                                    //  配合sizeof(变量)使用
{
    assert(el_sz > 0);  // 这个不用看

    zarray_t *za = (zarray_t*) calloc(1, sizeof(zarray_t)); // 创建zarray_t对象
    za->el_sz = el_sz; // 指定zarray中存储的元素大小
    return za; // 返回对象
}
zarray中的扩容规则,可以进行修改,减少扩容次数,加速算法
static inline void zarray_ensure_capacity(zarray_t *za, int capacity)
{
    assert(za != NULL); // 判断 不用看

    if (capacity <= za->alloc) // 判断 不用看
        return;

    while (za->alloc < capacity) {
        za->alloc *= 2; // 每次x2 进行扩容 与 C++ Java 类似
        if (za->alloc < 8)
            za->alloc = 8;
    }

    za->data = (char*) realloc(za->data, za->alloc * za->el_sz); // 重新分配
}
元素添加操作
static inline void zarray_add(zarray_t *za, const void *p) // p待添加的元素
{
    assert(za != NULL); // 不用看
    assert(p != NULL); // 不用看

    zarray_ensure_capacity(za, za->size + 1); // 确保够大

    memcpy(&za->data[za->size*za->el_sz], p, za->el_sz); // 将p指向的内存放入zarray后面 是一个copy操作 
    za->size++; // 元素数加1
}
元素获取
static inline void zarray_get(const zarray_t *za, int idx, void *p) // idx获取第几个元素,p存放获取的元素
{
    assert(za != NULL); // 不用看
    assert(p != NULL); // 不用看
    assert(idx >= 0); // 不用看
    assert(idx < za->size); // 不用看

    memcpy(p, &za->data[idx*za->el_sz], za->el_sz); // 拷贝到p中
}

图像结构体

typedef struct image_u8 image_u8_t; // 图像结构体
struct image_u8 // uint8_t 的图像 uchar Mat
{
    const int32_t width; // 图像的宽度
    const int32_t height; // 图像的高度
    const int32_t stride; // 图像的遍历方式 类似于Mat中的step 
    // 举个例子 1280x720 的图 width = 1280 
    //                       height = 720
    //                       stride = 1280
    // 使用时 image_u8_t* image;
    // 访问元素 image.buf[rowId*stride + colId]; // 访问(rowId, colId)元素
    uint8_t *buf; // 存储像素的空间
};

核心检测代码部分

最为核心的函数 **zarray_t *apriltag_detector_detect(apriltag_detector_t td, image_u8_t im_orig)

// apriltag_detector_t *td 这个是检测器的一些基本参数配置(你可以设置的参数)
// image_u8_t *im_orig 原始图像
zarray_t *apriltag_detector_detect(apriltag_detector_t *td, image_u8_t *im_orig)
{
	if (zarray_size(td->tag_families) == 0) { // 一个设置参数问题,如果没有对应的tag码就会报错
	    zarray_t *s = zarray_create(sizeof(apriltag_detection_t*));
	    debug_print("No tag families enabled\n");
	    return s;
    } // 基本上没用
	
	// 这个也是一个线程操作 根据用户设置的参数决定
	if (td->wp == NULL || td->nthreads != workerpool_get_nthreads(td->wp)) { 
	    workerpool_destroy(td->wp);
	    td->wp = workerpool_create(td->nthreads);
	    if (td->wp == NULL) {
	        // creating workerpool failed - return empty zarray
	        return zarray_create(sizeof(apriltag_detection_t*));
	    }
	}
	
	///
	// Step 1. 图像预处理环节
	image_u8_t *quad_im = im_orig;
	if (td->quad_decimate > 1) { // 如果使用者设置了这个参数,就会执行图像下采样(可以加速算法)
	    quad_im = image_u8_decimate(im_orig, td->quad_decimate); // 主要由这个函数实现,这个没什么困难的,可以自己看一下
	    // 作者在这里进行了一个改动,对下采样为1.5的情况进行了加权采样
	}

	// 如果设置了这个参数,就会进行高斯滤波操作,具体可以看图像处理的文章,这里不是关键
	if (td->quad_sigma != 0) { // 高斯滤波
    	...
    }
	
	if (td->debug) // 如果设置了debug就可以看处理后的结果图 这里需要转换 网上有pnm转化的工具,也可以自己写一下
    	image_u8_write_pnm(quad_im, "debug_preprocess.pnm");	
    // 最最重要的部分来了 (单独解释一下这个函数)
	zarray_t *quads = apriltag_quad_thresh(td, quad_im);
	
	... 后面就是解码部分,没怎么细看,感觉应该是将图像点通过Hinv转到原始三维平面点进行的快速解码操作
}

重点方法分析

*zarray_t quads = apriltag_quad_thresh(td, quad_im);

// apriltag_detector_t *td 检测器设置变量
// image_u8_t *im // 下采样/高斯后的图像
zarray_t *apriltag_quad_thresh(apriltag_detector_t *td, image_u8_t *im)
{
    
    // step 1. 图像"二值化"操作,其实按照作者的操作,一共保留了三个状态 0 255 127
    //  0: 黑
    //  255:白 
    //  127:默认不进行操作的值,这个值是由我们设定的黑白像素差来决定的,小于这个差就认为边界差距不大,不做处理(后面会讲)

    int w = im->width, h = im->height; // 这个就是保存图像的宽和高

    image_u8_t *threshim = threshold(td, im); // 这个就是黑白阈值操作(单独讲,见下面)
    int ts = threshim->stride;

    if (td->debug)
        image_u8_write_pnm(threshim, "debug_threshold.pnm");
    
    // step 2. 寻找连通域 学习一下connected_components 见下面
    unionfind_t* uf = connected_components(td, threshim, w, h, ts);

    // make segmentation image.
    if (td->debug) {
    	// 我删除了debug部分的代码,这是一个画图展示的部分,太占篇幅
    }
    zarray_t* clusters = gradient_clusters(td, threshim, w, h, ts, uf);

    if (td->debug) {
       	// 老规矩我删除了画图
    }

    image_u8_destroy(threshim);

    
    // step 3. 处理四边形部分我没怎么看,感觉优化不了啥 sorry

    zarray_t* quads = fit_quads(td, w, h, clusters, im);
    unionfind_destroy(uf);

    for (int i = 0; i < zarray_size(clusters); i++) {
        zarray_t *cluster;
        zarray_get(clusters, i, &cluster);
        zarray_destroy(cluster);
    }
    zarray_destroy(clusters);

    return quads;
}

threshold 黑白阈值操作

代码中部分assert操作我删除了,不然太长了。我基本保留了所有的代码。

// apriltag_detector_t *td 检测器配置变量
// image_u8_t *im 下采样/高斯滤波处理后的图像
image_u8_t *threshold(apriltag_detector_t *td, image_u8_t *im)
{
    int w = im->width, h = im->height, s = im->stride; // 取出图像的基本信息

    image_u8_t *threshim = image_u8_create_alignment(w, h, s); // 一个对齐操作,threshim 这个是最终的返回图
    const int tilesz = 4; // 这个是论文中所提到的4x4小块  也可以进行调参
    
    int tw = w / tilesz; // 计算宽上有多少个块
    int th = h / tilesz; // 计算高上有多少个块 总的块数为 tw * th,因为这里是进行int操作,可能会存在余数的情况,后面作者进行了操作 (请注意这些块不存在重叠)

    uint8_t *im_max = calloc(tw*th, sizeof(uint8_t)); // 申明总块数的操作 记录每一个块中的最大值
    uint8_t *im_min = calloc(tw*th, sizeof(uint8_t)); // 申明总块数的操作 记录每一个块中的最小值

    struct minmax_task *minmax_tasks = malloc(sizeof(struct minmax_task)*th); // 这是线程操作 不用太在意, 主要意思就是 申明了th行的任务数,进行操作
   	// 开始操作,按照行进行块的最大最小值搜索
    for (int ty = 0; ty < th; ty++) {
        minmax_tasks[ty].im = im; // 将相应的变量保存至对应的参数中
        minmax_tasks[ty].im_max = im_max;
        minmax_tasks[ty].im_min = im_min;
        minmax_tasks[ty].ty = ty;

        workerpool_add_task(td->wp, do_minmax_task, &minmax_tasks[ty]); 
        // 这个是线程任务,你主要关注第二个参数就行,这是个函数指针,do_minmax_task是一个函数,后面基本上所有的代码都是这种格式,熟悉它就行。
    }
    workerpool_run(td->wp); // 开始执行线程任务 也就是上面的循环加入的操作
    free(minmax_tasks);

   	// 对得到的tw * th个块的最大最小值矩阵进行3x3的最大最小滤波,目的就是为了处理上面谈到的4x4块之间没有重叠区域导致的边界问题(论文中有,可以大概看一下)
    if (1) {
        uint8_t *im_max_tmp = calloc(tw*th, sizeof(uint8_t)); // 保存更新后的最大最小值表
        uint8_t *im_min_tmp = calloc(tw*th, sizeof(uint8_t));

        struct blur_task *blur_tasks = malloc(sizeof(struct blur_task)*th);
        for (int ty = 0; ty < th; ty++) {
            blur_tasks[ty].im = im;
            blur_tasks[ty].im_max = im_max;
            blur_tasks[ty].im_min = im_min;
            blur_tasks[ty].im_max_tmp = im_max_tmp;
            blur_tasks[ty].im_min_tmp = im_min_tmp;
            blur_tasks[ty].ty = ty;

            workerpool_add_task(td->wp, do_blur_task, &blur_tasks[ty]); // 老规矩进行线程配置 ,主要看do_blur_task函数
        }
        workerpool_run(td->wp);
        free(blur_tasks);
        free(im_max);
        free(im_min);
        im_max = im_max_tmp;
        im_min = im_min_tmp;
    }
	
	// 二值化操作
    struct threshold_task *threshold_tasks = malloc(sizeof(struct threshold_task)*th);
    for (int ty = 0; ty < th; ty++) {
        threshold_tasks[ty].im = im; // 传入处理后的图
        threshold_tasks[ty].threshim = threshim; // 传入最终二值化的图
        threshold_tasks[ty].im_max = im_max; // 上一步计算的最大最小值矩阵
        threshold_tasks[ty].im_min = im_min;
        threshold_tasks[ty].ty = ty; // 第几行
        threshold_tasks[ty].td = td;

        workerpool_add_task(td->wp, do_threshold_task, &threshold_tasks[ty]); // 老规矩 看do_threshold_task
    }
    workerpool_run(td->wp);
    free(threshold_tasks);

    // 这里就是作者对不满足4x4大小的边缘部分的处理(这部分我没仔细看,我认为对整个代码的梳理不是很重要,如果谁看了可以给我解释一下,谢谢🤭)
    if (1) {
        for (int y = 0; y < h; y++) {
            int x0;
            if (y >= th*tilesz) {
                x0 = 0; 
            } else {
                x0 = tw*tilesz; 
            }

            // compute tile coordinates and clamp.
            int ty = y / tilesz;
            if (ty >= th)
                ty = th - 1;

            for (int x = x0; x < w; x++) {
                int tx = x / tilesz;
                if (tx >= tw)
                    tx = tw - 1;

                int max = im_max[ty*tw + tx];
                int min = im_min[ty*tw + tx];
                int thresh = min + (max - min) / 2;

                uint8_t v = im->buf[y*s+x];
                if (v > thresh)
                    threshim->buf[y*s+x] = 255;
                else
                    threshim->buf[y*s+x] = 0;
            }
        }
    }

    free(im_min);
    free(im_max);

	
    if (td->qtp.deglitch) {
    	.. 这部分代码按照作者的意思貌似没啥改善
    }

    return threshim; // 返回“二值化”图
}

至此我们介绍完了二值化函数

do_minmax_task方法
void do_minmax_task(void *p)
{
    const int tilesz = 4; // 4x4块
    struct minmax_task* task = (struct minmax_task*) p; // 强制转化类型
    int s = task->im->stride; 
    int ty = task->ty;
    int tw = task->im->width / tilesz;
    image_u8_t *im = task->im;

    for (int tx = 0; tx < tw; tx++) { // 由于每一线程执行一行搜索,所以这里只有按行进行的操作
        uint8_t max = 0, min = 255;

        for (int dy = 0; dy < tilesz; dy++) { // 这个就对4x4的小块进行的最大最小值查找 
            for (int dx = 0; dx < tilesz; dx++) {

                uint8_t v = im->buf[(ty*tilesz+dy)*s + tx*tilesz + dx];
                if (v < min)
                    min = v;
                if (v > max)
                    max = v;
            }
        }
        task->im_max[ty*tw+tx] = max; // 将当前的4x4的块中的最大值放入task变量中,传出函数
        task->im_min[ty*tw+tx] = min; // 将当前的4x4的块中的最小值放入task变量中,传出函数
    }
}
do_blur_task 我删除了一些基本操作(太长了)
void do_blur_task(void *p)
{
    const int tilesz = 4; 
    for (int tx = 0; tx < tw; tx++) { // 还是按行操作
        uint8_t max = 0, min = 255;

        for (int dy = -1; dy <= 1; dy++) { // -1 0 1 3个 
            if (ty+dy < 0 || ty+dy >= th) // 判断是否越界
                continue;
            for (int dx = -1; dx <= 1; dx++) { 
                if (tx+dx < 0 || tx+dx >= tw) // 判断是否越界
                    continue;

                uint8_t m = im_max[(ty+dy)*tw+tx+dx];
                if (m > max)
                    max = m;
                m = im_min[(ty+dy)*tw+tx+dx];
                if (m < min)
                    min = m;
            }
        }

        task->im_max_tmp[ty*tw + tx] = max; // 更新每一个块中的最大最小值
        task->im_min_tmp[ty*tw + tx] = min;
    }
}
do_threshold_task 这个比较重要没删减
void do_threshold_task(void *p)
{
    const int tilesz = 4;
	struct threshold_task* task = (struct threshold_task*) p;
	int ty = task->ty;
	int tw = task->im->width / tilesz;
	int s = task->im->stride;
	uint8_t *im_max = task->im_max;
	uint8_t *im_min = task->im_min;
	image_u8_t *im = task->im;
	image_u8_t *threshim = task->threshim;
	int min_white_black_diff = task->td->qtp.min_white_black_diff; // 这个参数很重要,是一个可以调整的参数,就是上面说的产生127的参数

    for (int tx = 0; tx < tw; tx++) { // 仍然只进行 行间操作
        int min = im_min[ty*tw + tx]; // 获取4x4块中的最大最小值 用于计算阈值
        int max = im_max[ty*tw + tx];
        
        if (max - min < min_white_black_diff) { // 如果(最大最小)像素之差小于参数黑白差,则将图像的这部分置为127
            for (int dy = 0; dy < tilesz; dy++) {
                int y = ty*tilesz + dy;

                for (int dx = 0; dx < tilesz; dx++) {
                    int x = tx*tilesz + dx;

                    threshim->buf[y*s+x] = 127;
                }
            }
            continue;
        }

        // 否则进行0 255二值化操作
        uint8_t thresh = min + (max - min) / 2; // 近似为 (max + min)/2

        for (int dy = 0; dy < tilesz; dy++) { // 4 x 4快中图像按阈值变为 0 或 255
            int y = ty*tilesz + dy;

            for (int dx = 0; dx < tilesz; dx++) {
                int x = tx*tilesz + dx;

                uint8_t v = im->buf[y*s+x];
                if (v > thresh) // 大于阈值变为 255
                    threshim->buf[y*s+x] = 255;
                else 
                    threshim->buf[y*s+x] = 0;
            }
        }
    }
}
connected_components 代码可能进行了部分删减
unionfind_t* connected_components(apriltag_detector_t *td, image_u8_t* threshim, int w, int h, int ts) {
    unionfind_t *uf = unionfind_create(w * h); // 并查集初始化 // 并查集的大小为 图像的大小(初始化)

    if (td->nthreads <= 1) { // 由于我只进行了单线程,我就这部分讲解,多线程情况没看 我删除了 sorry 
        do_unionfind_first_line(uf, threshim, w, ts); // 判断第一行的像素连通情况
        for (int y = 1; y < h; y++) {
            do_unionfind_line2(uf, threshim, w, ts, y); // 判断下面所有行的连通情况,这样的原因在于作者查找连通域时的设计方案
        }
    } else {
      	// 多线程配置
    }
    return uf;
}
unionfind_create 并查集初始化
// 并查集合并操作 很重要
// unionfind_t *uf 并查集  
// uint32_t aid 集合1
// uint32_t bid 集合2
static inline uint32_t unionfind_connect(unionfind_t *uf, uint32_t aid, uint32_t bid)
{
    uint32_t aroot = unionfind_get_representative(uf, aid); // 获取root
    uint32_t broot = unionfind_get_representative(uf, bid); 

    if (aroot == broot) // 如果两个集合为相同的集合就返回
        return aroot;

    uint32_t asize = uf->size[aroot] + 1; // 否则集合数目加1
    uint32_t bsize = uf->size[broot] + 1;

    if (asize > bsize) { // 向集合数目多的那个合并
        uf->parent[broot] = aroot;
        uf->size[aroot] += bsize;
        return aroot;
    } else {
        uf->parent[aroot] = broot;
        uf->size[broot] += asize;
        return broot;
    }
}

static inline uint32_t unionfind_get_representative(unionfind_t *uf, uint32_t id)
{
    uint32_t root = uf->parent[id]; 
    // 未初始化的节点
    if (root == 0xffffffff) { 
        uf->parent[id] = id; // 将对应的根设置为自己的Id
        return id;
    }

    // 找到自己的根
    while (uf->parent[root] != root) {
        root = uf->parent[root];
    }

    // 回溯修改所有遍历的上级节点的根为最终的根节点的Id,便于快速查找
    while (uf->parent[id] != root) {
        uint32_t tmp = uf->parent[id];
        uf->parent[id] = root;
        id = tmp;
    }

    return root;
}

static inline unionfind_t *unionfind_create(uint32_t maxid)
{
    unionfind_t *uf = (unionfind_t*) calloc(1, sizeof(unionfind_t)); // 并查集分配内存
    uf->maxid = maxid;
    uf->parent = (uint32_t *) malloc((maxid+1) * sizeof(uint32_t) * 2); // 分配内存给parent size 注意这里多分配了一个
    memset(uf->parent, 0xff, (maxid+1) * sizeof(uint32_t)); // 初始化所有的parent为0xff
    uf->size = uf->parent + (maxid+1); // size 指向了刚刚多分配的那个内存块
    memset(uf->size, 0, (maxid+1) * sizeof(uint32_t)); // 又给并查集中的size分配内存
    return uf;
}

do_unionfind_first_line

#define DO_UNIONFIND2(dx, dy) if (im->buf[(y + dy)*s + x + dx] == v) unionfind_connect(uf, y*w + x, (y + dy)*w + x + dx);

static void do_unionfind_first_line(unionfind_t *uf, image_u8_t *im, int w, int s)
{
    int y = 0;
    uint8_t v;

    for (int x = 1; x < w - 1; x++) { // 注意这里是从1开始的
        v = im->buf[y*s + x]; // 获取当前的元素 (0,1) (0,2) ,.....

        if (v == 127) // 如果元素为127 跳过
            continue;

        DO_UNIONFIND2(-1, 0); // 判断连接 (0,0)-(0,1) (0,1)-(0,2), ........
        					  // (im->buf[(y + dy)*s + x + dx] == v) 如果相等就将两个集合连接起来
    }
}

do_uninfin_line2 (我删除了assert部分代码)

static void do_unionfind_line2(unionfind_t *uf, image_u8_t *im, int w, int s, int y)
{
    uint8_t v_m1_m1;
    uint8_t v_0_m1 = im->buf[(y - 1)*s];
    uint8_t v_1_m1 = im->buf[(y - 1)*s + 1];
    uint8_t v_m1_0;
    uint8_t v = im->buf[y*s];

    for (int x = 1; x < w - 1; x++) { // 其实这部分看懂只要在纸上画画就知道了,我没办法讲的太详细,但是这里面其实是有点分类不全的。
        v_m1_m1 = v_0_m1;
        v_0_m1 = v_1_m1;
        v_1_m1 = im->buf[(y - 1)*s + x + 1];
        v_m1_0 = v;
        v = im->buf[y*s + x];

        if (v == 127)
            continue;

        // (-1, -1)    (0, -1)    (1, -1)
        // (-1, 0)    
        DO_UNIONFIND2(-1, 0);

        if (x == 1 || !((v_m1_0 == v_m1_m1) && (v_m1_m1 == v_0_m1))) {
            DO_UNIONFIND2(0, -1);
        }

        if (v == 255) {
            if (x == 1 || !(v_m1_0 == v_m1_m1 || v_0_m1 == v_m1_m1) ) {
                DO_UNIONFIND2(-1, -1);
            }
            if (!(v_0_m1 == v_1_m1)) {
                DO_UNIONFIND2(1, -1); // 具体看下图
            }
        }
    }
}

至此连通域解释完毕
连通区域图
连通域Bug

gradient_clusters
// 这是一个会用到的类
struct cluster_hash
{
    uint32_t hash; // hash值 hash值是不唯一
    uint64_t id; // id值 需要注意 id值是唯一的
    zarray_t* data; // 这个data具体是指struct pt
};

zarray_t* gradient_clusters(apriltag_detector_t *td, image_u8_t* threshim, int w, int h, int ts, unionfind_t* uf) {
    zarray_t* clusters; // 初始化一个聚类的zarray_t对象 还未创建
    int nclustermap = 0.2*w*h;

    int sz = h - 1;
    int chunksize = 1 + sz / (APRILTAG_TASKS_PER_THREAD_TARGET * td->nthreads); // 这里的参数APRILTAG_TASKS_PER_THREAD_TARGET也是可以调整的,主要用于调整每次的任务数量,任务越多内存可能就越散。
    struct cluster_task *tasks = malloc(sizeof(struct cluster_task)*(sz / chunksize + 1));
    int ntasks = 0;

    for (int i = 1; i < sz; i += chunksize) {
        // 每一个任务执行[y0, y1). 
        tasks[ntasks].y0 = i;
        tasks[ntasks].y1 = imin(sz, i + chunksize);
        tasks[ntasks].w = w;
        tasks[ntasks].s = ts;
        tasks[ntasks].uf = uf;
        tasks[ntasks].im = threshim;
        tasks[ntasks].nclustermap = nclustermap/(sz / chunksize + 1);
        tasks[ntasks].clusters = zarray_create(sizeof(struct cluster_hash*)); // 创建一个cluster_hash的vector对象

        workerpool_add_task(td->wp, do_cluster_task, &tasks[ntasks]); // 老规矩 do_cluster_task 这个很抽象 可以想象成将连通域的相邻边界记录下来
        ntasks++;
    }

    workerpool_run(td->wp);
	
	// 下面的操作是将上面多个tasks中的分开(每个任务开了一个独立的hash表)的cluster_hash进行合并操作,能够合并的原因在于使用了同样的全局唯一ID
    zarray_t** clusters_list = malloc(sizeof(zarray_t *)*ntasks); // 申明ntasks个zarray_t对象, 这个ntasks是前面记录你开了几个任务进行的,减少这个可以加速算法
    for (int i = 0; i < ntasks; i++) {
        clusters_list[i] = tasks[i].clusters; // 将do_cluster_task中获得的hash表换个地方存
    }

    int length = ntasks;
    // 此处也是很有技巧性的操作,有归并排序的感觉 不断除2 合并
    while (length > 1) { 
        int write = 0;
        for (int i = 0; i < length - 1; i += 2) { 
            clusters_list[write] = merge_clusters(clusters_list[i], clusters_list[i + 1]); // 合并完成了可以进行覆盖
            write++; 
        }

        if (length % 2) { // 当长度为奇数时, 进行如下的操作, 单独加入
            clusters_list[write] = clusters_list[length - 1];
        }

        length = (length >> 1) + length % 2; // 除2 加 余数
    }
	// 完成后所有的数据全部集中到了clusters_list[0]中
	
    clusters = zarray_create(sizeof(zarray_t*));
    zarray_ensure_capacity(clusters, zarray_size(clusters_list[0])); // 将元素为struct cluster_hash*的zarray_t*放入其中
    for (int i = 0; i < zarray_size(clusters_list[0]); i++) {
        struct cluster_hash** hash;
        zarray_get_volatile(clusters_list[0], i, &hash);
        zarray_add(clusters, &(*hash)->data); // 只保留hash中的data数据,也就是边的数据 C真的太强大了,全靠的zarray_t*指针,data也时zarray_t*类型
        free(*hash);
    }
    zarray_destroy(clusters_list[0]);
    free(clusters_list);
    free(tasks);
    return clusters;
}

static void do_cluster_task(void *p)
{
    struct cluster_task *task = (struct cluster_task*) p;
    do_gradient_clusters(task->im, task->s, task->y0, task->y1, task->w, task->nclustermap, task->uf, task->clusters); // 主要执行这个函数
}

struct uint64_zarray_entry // 这是一个链表结构 这其实是一个hash表
{
    uint64_t id;
    zarray_t *cluster;

    struct uint64_zarray_entry *next;
};

zarray_t* do_gradient_clusters(image_u8_t* threshim, int ts, int y0, int y1, int w, int nclustermap, unionfind_t* uf, zarray_t* clusters) {
    struct uint64_zarray_entry **clustermap = calloc(nclustermap, sizeof(struct uint64_zarray_entry*)); // 此处申明了一个二维指针 nclustermap个

    int mem_chunk_size = 2048;
    struct uint64_zarray_entry** mem_pools = malloc(sizeof(struct uint64_zarray_entry *)*(1 + 2 * nclustermap / mem_chunk_size)); // SmodeTech: avoid memory corruption when nclustermap < mem_chunk_size 这个不是很清楚 sorry 并不是非常影响代码阅读 有点类似先开辟了这块空间待使用
    int mem_pool_idx = 0;
    int mem_pool_loc = 0;
    mem_pools[mem_pool_idx] = calloc(mem_chunk_size, sizeof(struct uint64_zarray_entry));

    for (int y = y0; y < y1; y++) { // 处理[y0, y1)的任务
        bool connected_last = false;
        for (int x = 1; x < w-1; x++) { // 开始遍历

            uint8_t v0 = threshim->buf[y*ts + x]; // 当前位置的像素 0 127 255
            if (v0 == 127) { // 如果是127跳过 
                connected_last = false;
                continue;
            }

            uint64_t rep0 = unionfind_get_representative(uf, y*w + x); // 使用并查集找到当前位置的根节点,也就是所属集合
            if (unionfind_get_set_size(uf, rep0) < 25) { // 如果所属集合中的个数少于一定的像素,则最直接跳过, 这个25也是一个可调整的参数
                connected_last = false;
                continue;
            }
            bool connected; // 这个默认为false
            // 下面是一个宏定义的函数 也是最核心的一个函数 主要目的是获得两个集合的边界的位置与梯度
#define DO_CONN(dx, dy)                                                 \
            if (1) {                                                    \
            	/*取出当前像素邻近的某一个像素点*/                        \                    
                uint8_t v1 = threshim->buf[(y + dy)*ts + x + dx];       \ 
                														\
                /*判断这两个相邻像素的值是否相差255,也就是存在梯度*/      \
                if (v0 + v1 == 255) {                                   \
                	/*找到它的根节点*/
                    uint64_t rep1 = unionfind_get_representative(uf, (y + dy)*w + x + dx); \
                    /*这个也是可调整的参数,如果数量太少就舍弃*/
                    if (unionfind_get_set_size(uf, rep1) > 24) {        \
                        uint64_t clusterid;                                 \
                        /*将两个集合的Id/根 进行从大到小的组合, 具体来说是将大的数放在前32位,小大数放在后32位, 注意这个clusterid是唯一,只取决于相邻的两个连通域id*/
                        if (rep0 < rep1)                                    \
                            clusterid = (rep1 << 32) + rep0;                \
                        else                                                \
                            clusterid = (rep0 << 32) + rep1;                \
                                                                            \
                        /*使用lousy hash 函数,按这个规则放入某一个数组中(这步求的就是数组下标)*/                       \
                        uint32_t clustermap_bucket = u64hash_2(clusterid) % nclustermap; \
                        /*使用前面申明的二维指针,此时该对象还未初始化,下面就是初始化操作以及添加操作*/
                        struct uint64_zarray_entry *entry = clustermap[clustermap_bucket]; \
                        /*这是一个链表操作,判断hash表中是否存在当前的id,终止条件为entry为空或者存在相同的id*/
                        while (entry && entry->id != clusterid) {           \
                            entry = entry->next;                            \
                        }                                                   \
                                                                            \
                        /*这个地方只有在为空时才会进入*/
                        /*这部分主要是初始化*/
                        if (!entry) {                                       \
                            if (mem_pool_loc == mem_chunk_size) {           \
                            /*在之前申明的2048个块中周期性取*/
                                mem_pool_loc = 0;                           \
                                mem_pool_idx++;                             \
                                mem_pools[mem_pool_idx] = calloc(mem_chunk_size, sizeof(struct uint64_zarray_entry)); \
                            }                                               \
                            entry = mem_pools[mem_pool_idx] + mem_pool_loc; \
                            mem_pool_loc++;                                 \
                                                                            \
                            /*这个标识是唯一的*/
                            entry->id = clusterid;                          \
                            /*这行代码是说entry中的zarray_t数组是存储struct pt类型的,这个类型内部有5个变量, 当前的x,y像素位置, 当前的位置处的gx gy梯度, 以及斜率*/
                            /*也就说凡是相同的id的边全部集中到了这个entry中*/
                            entry->cluster = zarray_create(sizeof(struct pt)); \
                            /*下面两行是说每次都将最新的entry作为链表的头,可以认为是减少查找次数,因为相邻的大多数像素都是在一个集合中的*/
                            entry->next = clustermap[clustermap_bucket];    \
                            clustermap[clustermap_bucket] = entry;          \
                        }                                                   \
                                                                            \
                        /*注意这里的点是2x 2y实际应该是x y后面处理会进行处理*/
                        struct pt p = { .x = 2*x + dx, .y = 2*y + dy, .gx = dx*((int) v1-v0), .gy = dy*((int) v1-v0)}; \
                        /*将信息添加进聚类的边中*/
                        zarray_add(entry->cluster, &p);                     \
                        connected = true;                                   \
                    }                                                   \
                }                                                       \
            }
            DO_CONN(1, 0); // 开始查找边 右侧像素
            DO_CONN(0, 1);  // 下方像素

            if (!connected_last) { // 减少遍历次数, 可以自己画图看看
                DO_CONN(-1, 1);
            }
            connected = false;
            DO_CONN(1, 1);
            connected_last = connected;
        }
    }
#undef DO_CONN

    for (int i = 0; i < nclustermap; i++) { // 遍历上面保存的所有边的数据
        int start = zarray_size(clusters);
        for (struct uint64_zarray_entry *entry = clustermap[i]; entry; entry = entry->next) // 遍历hash表中的某一个链表
        {
            struct cluster_hash* cluster_hash = malloc(sizeof(struct cluster_hash)); 
            // 创建cluster_hash这个类用于保存上面获取的边的信息,同时放入cluster中,传出函数,以便后续使用
            cluster_hash->hash = u64hash_2(entry->id) % nclustermap; // 获取hash 实际上就是i
            cluster_hash->id = entry->id; // 这个id根据上面的解释是唯一的
            cluster_hash->data = entry->cluster; // 这个data是上面的entry->cluster = zarray_create(sizeof(struct pt));
            zarray_add(clusters, &cluster_hash); // 将信息加入后传出函数
        }
        int end = zarray_size(clusters);

        // 对id值进行冒泡排序,可以尝试其他排序方式,不知道怎么样,没试过🙂
        int n = end - start;
        for (int j = 0; j < n - 1; j++) {
            for (int k = 0; k < n - j - 1; k++) {
                struct cluster_hash** hash1;
                struct cluster_hash** hash2;
                zarray_get_volatile(clusters, start + k, &hash1);
                zarray_get_volatile(clusters, start + k + 1, &hash2);
                if ((*hash1)->id > (*hash2)->id) {
                    struct cluster_hash tmp = **hash2;
                    **hash2 = **hash1;
                    **hash1 = tmp;
                }
            }
        }
    }
    for (int i = 0; i <= mem_pool_idx; i++) {
        free(mem_pools[i]);
    }
    free(mem_pools);
    free(clustermap);

    return clusters;
}

至此完成了梯度聚类部分操作

merge_clusters 合并操作
// zarray_t* c1 具体指向 struct cluster_hash*
// zarray_t* c2 具体指向 struct cluster_hash*
// 合并有三种情况 c1 c2 存在交集 进行合并
// 两个完全没有交集
zarray_t* merge_clusters(zarray_t* c1, zarray_t* c2) {
    zarray_t* ret = zarray_create(sizeof(struct cluster_hash*)); // 返回值,存储合并的结果
    zarray_ensure_capacity(ret, zarray_size(c1) + zarray_size(c2));

    int i1 = 0;
    int i2 = 0;
    int l1 = zarray_size(c1);
    int l2 = zarray_size(c2);

    while (i1 < l1 && i2 < l2) { // 只有这两个都存在值的时候才进行操作
        struct cluster_hash** h1;
        struct cluster_hash** h2;
        zarray_get_volatile(c1, i1, &h1);
        zarray_get_volatile(c2, i2, &h2);

		// hash值相同id不一定相同, 但是id相同hash一定相同
        if ((*h1)->hash == (*h2)->hash && (*h1)->id == (*h2)->id) {
            zarray_add_range((*h1)->data, (*h2)->data, 0, zarray_size((*h2)->data)); // cluster_hash中的data存的是struct pt
            zarray_add(ret, h1);
            i1++;
            i2++;
            zarray_destroy((*h2)->data);
            free(*h2);
        } else if ((*h2)->hash < (*h1)->hash || ((*h2)->hash == (*h1)->hash && (*h2)->id < (*h1)->id)) { // 这边可以这样操作的原因是什么? 还记得我们在梯度聚类中的冒泡排序吗?
            zarray_add(ret, h2);
            i2++;
        } else {
            zarray_add(ret, h1);
            i1++;
        }
    }

    zarray_add_range(ret, c1, i1, l1); // 将剩余的对象全部加入ret中
    zarray_add_range(ret, c2, i2, l2);

    zarray_destroy(c1);
    zarray_destroy(c2);

    return ret;
}

合并操作示意图
hash表合并操作

至此完成了梯度聚类所有的代码讲解

后面的部分是判断是否符合tag码的设计要求,以及解码操作没怎么看了

看完源码后,我只能想说C是真的美!

  • 4
    点赞
  • 7
    收藏
    觉得还不错? 一键收藏
  • 1
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值