【代码阅读】darknet源码阅读（四）:matrix.h 和 matrix.c

参考文献依然是放前面：https://blog.youkuaiyun.com/caicaiatnbu/category_9096319.html

darknet版本： https://github.com/AlexeyAB/darknet，与原始的版本还是有一点区别的。

进入代码：

1.matrix.h

#ifndef MATRIX_H
#define MATRIX_H
#include "darknet.h"

//typedef struct matrix{
//    int rows, cols;
//    float **vals;
//} matrix;
//与原版不同，这里添加定义了结构体model，用到了matrix啊，这个结构体定义在darknet.h里，结构如上注释
typedef struct {
    int *assignments;
    matrix centers;
} model;

#ifdef __cplusplus
extern "C" {
#endif
//函数在matrix.c里看
model do_kmeans(matrix data, int k);
matrix make_matrix(int rows, int cols);
void free_matrix(matrix m);
void print_matrix(matrix m);

matrix csv_to_matrix(char *filename);
void matrix_to_csv(matrix m);
matrix hold_out_matrix(matrix *m, int n);
float matrix_topk_accuracy(matrix truth, matrix guess, int k);
void matrix_add_matrix(matrix from, matrix to);
void scale_matrix(matrix m, float scale);
matrix resize_matrix(matrix m, int size);

float *pop_column(matrix *m, int c);

#ifdef __cplusplus
}
#endif
#endif

2.matrix.c

#include "matrix.h"
#include "utils.h"
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <assert.h>
#include <math.h>


//释放矩阵m的空间
void free_matrix(matrix m)
{
    int i;
    for(i = 0; i < m.rows; ++i) free(m.vals[i]);//循环逐行释放空间
    free(m.vals);
}

//具体是什么操作，要看使用位置，等我遇见了再来看
//不清楚truth和guess具体指代值
float matrix_topk_accuracy(matrix truth, matrix guess, int k)
{
    int* indexes = (int*)xcalloc(k, sizeof(int));//临时空间
    int n = truth.cols;//truth矩阵的列
    int i,j;
    int correct = 0;
    for(i = 0; i < truth.rows; ++i){
        //找到guess矩阵中，大小排在前k个的值，并返回前k个值对应的index放在indexes中
        top_k(guess.vals[i], n, k, indexes);//top_k定义在utils.c中
        //遍历top_k的数据
        for(j = 0; j < k; ++j){
            int class_id = indexes[j];
            if(truth.vals[i][class_id]){//对应 truth.vals[i][class] 位置是否非0
                ++correct;
                break;
            }
        }
    }
    free(indexes);// 释放index
    return (float)correct/truth.rows;// 返回比例。
}

//矩阵与常数乘法操作,将矩阵m中每个元素都放大scale倍
void scale_matrix(matrix m, float scale)
{
    int i,j;
    for(i = 0; i < m.rows; ++i){
        for(j = 0; j < m.cols; ++j){
            m.vals[i][j] *= scale;
        }
    }
}

//矩阵的行数和列数进行resize操作,resize矩阵大小是size * size
matrix resize_matrix(matrix m, int size)
{
    int i;
    if (m.rows == size) return m;//行数恰好等于size大小时，不做缩放
    if (m.rows < size) {//行数小于size时
        //xrealloc为重构原有存储空间函数，定义在utils.c中
        /*
        函数功能：动态调整内存，先判断当前的指针是否有足够的连续空间，如果有，则扩大mem_address指向的地址，
        并且将mem_address返回，如果空间不够，先按照newsize指定的大小分配空间，将原有数据全部拷贝到新分配的内存区域，
        而后对原来meme_address所指向的内存区域进行释放【这里是自动释放，不需要手动释放】，同时返回新分配内存区域的首地址。
        如果失败则返回空指针NULL；
*/
        m.vals = (float**)xrealloc(m.vals, size * sizeof(float*)); // 重新申请内存空间
        for (i = m.rows; i < size; ++i) {
            // 每一行的存储空间也要重新申请，每一行保存size 个数据
            m.vals[i] = (float*)xcalloc(m.cols, sizeof(float));
        }
    } else if (m.rows > size) {// 调整后矩阵行数减少
        for (i = size; i < m.rows; ++i) {// 释放多余的存储空间
            free(m.vals[i]);
        }
        m.vals = (float**)xrealloc(m.vals, size * sizeof(float*));// 分配每一行存储空间，每一行保存size 个数据
    
    }
    m.rows = size;
    return m; // m是一个size*size 大小的矩阵
}

//矩阵加法
void matrix_add_matrix(matrix from, matrix to)
{
    assert(from.rows == to.rows && from.cols == to.cols);
    int i,j;
    for(i = 0; i < from.rows; ++i){
        for(j = 0; j < from.cols; ++j){
            to.vals[i][j] += from.vals[i][j];
        }
    }
}


//创建大小一定的矩阵
matrix make_matrix(int rows, int cols)
{
    int i;
    matrix m;
    m.rows = rows;
    m.cols = cols;
    m.vals = (float**)xcalloc(m.rows, sizeof(float*));//创建空间
    for(i = 0; i < m.rows; ++i){
        m.vals[i] = (float*)xcalloc(m.cols, sizeof(float));//创建空间
    }
    return m;
}

//从矩阵中采样m行数据，返回采样后的结果；
matrix hold_out_matrix(matrix *m, int n)
{
    int i;
    matrix h;
    h.rows = n;
    h.cols = m->cols;
    h.vals = (float**)xcalloc(h.rows, sizeof(float*));//创建返回的采样矩阵空间
    for(i = 0; i < n; ++i){//循环需要的采样次数
        int index = rand()%m->rows;//随机采样一行
        h.vals[i] = m->vals[index];//赋值
        m->vals[index] = m->vals[--(m->rows)];// 把最后一行的数据覆盖到 index行上，填充？
    
    }
    return h;
}

//获取矩阵中某一列数据，并把该列删除掉，同时返回删除的值
float *pop_column(matrix *m, int c)
{
    float* col = (float*)xcalloc(m->rows, sizeof(float));//创建返回数组
    int i, j;
    for(i = 0; i < m->rows; ++i){//按行循环，删除该行该列的值
        col[i] = m->vals[i][c];//存储要被删掉的值
        for(j = c; j < m->cols-1; ++j){
            m->vals[i][j] = m->vals[i][j+1];//将删除掉的列后面的列往前移动，补满
        }
    }
    --m->cols;//将列数－1
    return col;//返回删除值
}

//读文件，加载矩阵
matrix csv_to_matrix(char *filename)
{
    FILE *fp = fopen(filename, "r");
    if(!fp) file_error(filename);//判断文件是否能读

    matrix m;//返回的矩阵
    m.cols = -1;

    char *line;

    int n = 0;
    int size = 1024;//初步设定矩阵有1024行
    m.vals = (float**)xcalloc(size, sizeof(float*));//分配空间
    while((line = fgetl(fp))){//按行读取字符数据
        //fgetl和count_fields定义在uitls.c中
        //count_fields：统计取字符串中有多少空格和',', 遇到第一个'\0'字符结束
        //','和空格作为字符的分割，也就是统计存储时矩阵需要多少列
        if(m.cols == -1) m.cols = count_fields(line);
        if(n == size){ //超过设定行数，则扩充矩阵存储空间
            size *= 2;
            m.vals = (float**)xrealloc(m.vals, size * sizeof(float*));
        }
        // parse_fields定义在uitls.c中
        //解析字符数组中m.cols个float小数,返回是一个float类型指针
        m.vals[n] = parse_fields(line, m.cols);
        free(line);
        ++n;
    }
    // 之前开辟的空间是1024的整数倍，此时需要根据n来实际分配内存空间。
    m.vals = (float**)xrealloc(m.vals, n * sizeof(float*));
    m.rows = n;
    return m;// 返回结果
}
//可视化矩阵m,打印
void matrix_to_csv(matrix m)
{
    int i, j;

    for(i = 0; i < m.rows; ++i){
        for(j = 0; j < m.cols; ++j){
            if(j > 0) printf(",");
            printf("%.17g", m.vals[i][j]);//自动选择合适的表示法输出
        }
        printf("\n");
    }
}

//可视化打印矩阵m
void print_matrix(matrix m)
{
    int i, j;
    printf("%d X %d Matrix:\n",m.rows, m.cols);//打印行和列数
    printf(" __");
    for(j = 0; j < 16*m.cols-1; ++j) printf(" ");
    printf("__ \n");

    printf("|  ");
    for(j = 0; j < 16*m.cols-1; ++j) printf(" ");
    printf("  |\n");

    for(i = 0; i < m.rows; ++i){
        printf("|  ");
        for(j = 0; j < m.cols; ++j){
            printf("%15.7f ", m.vals[i][j]);//打印浮点数内容
        }
        printf(" |\n");
    }
    printf("|__");
    for(j = 0; j < 16*m.cols-1; ++j) printf(" ");
    printf("__|\n");
}

matrix make_matrix(int rows, int cols);
//先申明，防止找不到函数
void copy(float *x, float *y, int n);
float dist(float *x, float *y, int n);
int *sample(int n);

//找到离datum最近的中心点，返回其位置
int closest_center(float *datum, matrix centers)
{
    int j;
    int best = 0;
    float best_dist = dist(datum, centers.vals[best], centers.cols);
    for (j = 0; j < centers.rows; ++j) {
        float new_dist = dist(datum, centers.vals[j], centers.cols);
        if (new_dist < best_dist) {
            best_dist = new_dist;//值越小就交换最优值
            best = j;
        }
    }
    return best;
}

//closest_center函数找到离datum最近的中心点位置，返回他们之间的距离
float dist_to_closest_center(float *datum, matrix centers)
{
    int ci = closest_center(datum, centers);
    return dist(datum, centers.vals[ci], centers.cols);
}

//kmeans聚类，循环样本，判断data中的每一个样本距离哪个centers更近，就属于哪一个centers的类别，并存放到assignments中去
//先看最后的函数do_kmeans,调用的这个函数
int kmeans_expectation(matrix data, int *assignments, matrix centers)
{
    int i;
    int converged = 1;//判断是否停止更新的标识符
    for (i = 0; i < data.rows; ++i) {//循环每一组样本的宽高
        //data中的宽高与初始的centers宽高做距离计算，获得最近的距离的index，
        //也就是获得对应样本的类别（距离哪个center更近，就属于哪个center的类别范围中）
        int closest = closest_center(data.vals[i], centers);
        //如果获得的类别，与原有存储中的类别不一样，则converged =0，更新标识符
        if (closest != assignments[i]) converged = 0;
        
        assignments[i] = closest;//更新相应的assignments（类别）值
    }
    return converged;返回更新的判断标识符
}


//看起来这里是kmeans生成anchor，
//assignments为m*1的数组，m为data的行数，存放的是data样本中每一个样本的类别
//data需要用到的训练样本的宽高
//centers初始聚类的宽高
void kmeans_maximization(matrix data, int *assignments, matrix centers)
{
    //初始化矩阵，为备份做准备
    matrix old_centers = make_matrix(centers.rows, centers.cols);
    int i, j;
    int *counts = (int*)xcalloc(centers.rows, sizeof(int));//分配空间

    //第一个循环，备份原有的center
    for (i = 0; i < centers.rows; ++i) {
        for (j = 0; j < centers.cols; ++j) {
            old_centers.vals[i][j] = centers.vals[i][j];//将center中的值放到old中
            centers.vals[i][j] = 0;//center清空
        }
    }
    
    //第二个循环，计算每个类别包含样本的样本值总和，和每个类别中包含的总样本数
    //counts中存放的是每一个类别包含的样本总数
    //centers这个循环后，存放的是该center类别包含的样本宽高值的总和
    for (i = 0; i < data.rows; ++i) {//按行循环所有的样本
        ++counts[assignments[i]];//
        for (j = 0; j < data.cols; ++j) {//
            centers.vals[assignments[i]][j] += data.vals[i][j];
        }
    }

    //第三个循环，计算新的center值=样本宽高总值/类别包含样本总数
    for (i = 0; i < centers.rows; ++i) {
        if (counts[i]) {
            for (j = 0; j < centers.cols; ++j) {
                //对每个类别刚刚计算的值除总数，计算平均值，更新centers
                centers.vals[i][j] /= counts[i];
            }
        }
    }

    //第四个循环，补上这一次没有更新的center值，将备份的值再赋给新center
    for (i = 0; i < centers.rows; ++i) {
        for (j = 0; j < centers.cols; ++j) {
            //如果有没计算的值，把原来的值再放回去
            if(centers.vals[i][j] == 0) centers.vals[i][j] = old_centers.vals[i][j];
        }
    }
    free(counts);
    free_matrix(old_centers);//释放临时空间
}


//随机生成centers
void random_centers(matrix data, matrix centers) {
    int i;
    int *s = sample(data.rows);//获得打乱的样本行数index值
    for (i = 0; i < centers.rows; ++i) {//循环centers矩阵
        copy(data.vals[s[i]], centers.vals[i], data.cols);//将按照打乱的index将data中的值放入centers中，获得初始随机的anchor值
    }
    free(s);
}


//利用sample生成随机的数组，数组的内容是数组对应的index
int *sample(int n)
{
    int i;
    int* s = (int*)xcalloc(n, sizeof(int));//构建返回数组空间
    for (i = 0; i < n; ++i) s[i] = i;//赋初值,值为数组index
    for (i = n - 1; i >= 0; --i) {
        int swap = s[i];
        int index = rand() % (i + 1);
        s[i] = s[index];
        s[index] = swap;//随机交换index,打乱存放的index
    }
    return s;//返回
}

//计算两个box之间iou距离,传入的x,y分别对应box的宽高
float dist(float *x, float *y, int n)
{
    //printf(" x0 = %f, x1 = %f, y0 = %f, y1 = %f \n", x[0], x[1], y[0], y[1]);
    float mw = (x[0] < y[0]) ? x[0] : y[0];
    float mh = (x[1] < y[1]) ? x[1] : y[1];
    float inter = mw*mh;
    float sum = x[0] * x[1] + y[0] * y[1];
    float un = sum - inter;
    float iou = inter / un;
    return 1 - iou;
}

//赋值数组值
void copy(float *x, float *y, int n)
{
    int i;
    for (i = 0; i < n; ++i) y[i] = x[i];
}

//做kmeans聚类
//data应该是m*2的矩阵，内容是训练样本box的宽，高
//k是需要的anchor数量，一般6或9
model do_kmeans(matrix data, int k)
{
    matrix centers = make_matrix(k, data.cols);//初始化矩阵

    //assignments将存放的是data样本中每一个样本的类别
    //（就是距离哪个centers更近，就属于哪个类别）
    int* assignments = (int*)xcalloc(data.rows, sizeof(int));//构建数组空间，大小为m*sizeof(int),m是我假设的data的行数，

    //smart_centers(data, centers);
    //获得随机的center(即anchor)值
    random_centers(data, centers);  // IoU = 67.31% after kmeans

    /*
    // IoU = 63.29%, anchors = 10,13,  16,30,  33,23,  30,61,  62,45,  59,119,  116,90,  156,198,  373,326
    centers.vals[0][0] = 10; centers.vals[0][1] = 13;
    centers.vals[1][0] = 16; centers.vals[1][1] = 30;
    centers.vals[2][0] = 33; centers.vals[2][1] = 23;
    centers.vals[3][0] = 30; centers.vals[3][1] = 61;
    centers.vals[4][0] = 62; centers.vals[4][1] = 45;
    centers.vals[5][0] = 59; centers.vals[5][1] = 119;
    centers.vals[6][0] = 116; centers.vals[6][1] = 90;
    centers.vals[7][0] = 156; centers.vals[7][1] = 198;
    centers.vals[8][0] = 373; centers.vals[8][1] = 326;
    */

    // range centers [min - max] using exp graph or Pyth example
    //如果anchor数量为1，则通过kmeans_maximization函数获得anchor值
    //不需要使用kmeans_expectation函数先给data中的样本分类，因为只有一个类别。
    if (k == 1) kmeans_maximization(data, assignments, centers);
    int i;
    
    //如果k>1哦，就需要使用kmeans_expectation函数先给data中的样本分类，将类别存放到assignments中
    //如果kmeans_expectation返回0，也就是说最近邻值有更新，就继续循环，否则停止更新。
    for(i = 0; i < 1000 && !kmeans_expectation(data, assignments, centers); ++i) {
        kmeans_maximization(data, assignments, centers);
    }
    printf("\n iterations = %d \n", i);
    model m;
    m.assignments = assignments;
    m.centers = centers;
    return m;
}