[人工智能]5万字用纯C语言从零开始实现人脸检测

前言

读者你好！ 这是我第一次写博客，也是我第一次对图像处理和人脸检测进行研究，所以这篇文章会有许多错误，请读者多多包含。但正如标题所述，是用C从读取二进制图像文件开始，到基本的图像处理，再到实现人脸检测（传统VJ算法），过程中并不需要下载任何第三方库。这篇文章记录了作者一整个的学习流程。且文中会毫无保留的展示所有功能的实现，和其相应的代码。
注意：1.每个知识点间的代码都是相互独立的，可直接复制进行测试运行。

一. BMP文件的处理

1.1 为什么选择BMP文件

BMP文件是未经压缩的文件，它保存了原图的所有的像素信息，文件体积最大，文件所包含的信息也最多，BMP的文件结构相对与JPG和PNG简单很多，易于处理和学习。（网络上一般以JPG格式的图片为主，所以使用时要进行格式转换，建议用电脑自带的画图工具转换。有些格式转换工具可能会加入一些其它的字节信息，破坏其原本结构）

1.2 BMP文件结构的介绍

1.）bmp文件的组成
bmp文件 = bmp头部（bmp文件头（14字节）+位图信息头（40字节））+位图数据 （bmp文件中的数据一个字节为一组，正向写入倒向读出）

2.）bmp文件头（14字节）
a.(2)424D---->‘BM’,表示这是Windows支持的位图格式。
b.(4)B6EA0500---->0050EAB6---->整个文件的大小（正向写入倒向读出）
c.(2)0000---->保留段，通常为0
d.(2)0000---->保留段，通常为0
e.(4)36000000---->从文件头到位图数据需偏移的字节数，即整个bmp头部的大小（54字节）

2.）位图信息头（40字节）
a.(4)28000000---->00000028---->位图信息头的大小（40）
b.(4)E9000000---->000000E9---->图片宽度（正向写入倒向读出）
c.(4)A0010000---->000001A0 ---->图片高度（正向写入倒向读出）
d.(2)0100->0001---->颜色平面数，即调色盘数，恒等于1
d.(2)2000->0020---->图片颜色的位数（32）
f.(4)00000000---->图像压缩比，不压缩则为
g.(4)00000000---->图像的大小，因为使用BI_RGB，所以设置为0
h.(4)00000000---->水平分辨率，缺省为0
i.(4)00000000---->垂直分辨率，缺省为0
j.(4)00000000---->说明本位图实际使用调色盘的颜色索引数，0表示全部
k.(4)00000000---->说明本位图重要调色盘的颜色索引数，0表示全都重要

1.3 BMP文件的读取和保存

代码说明与注意事项：
1. 代码中有详细的注释，核心就是C语言对二进制文件的读取和写入，但由于结构体对齐问题的存在，所以我选择最保险的，一条一条读取和写入的方式。
2. 虽然目前只能读取BMP文件，但还是保留了JPG的接口空位，后续可对其功能进行拓展。
3. BMP有32位24位16位8位，该程序只对前两个进行了处理
4. 计算机读取能被4整除大小的二进制文件时速度较快，所以读取24位图像的二进制文件时，每一行都会有补零操作。（每个图像每行补零的字节个数是固定的，且最多补3个字节）

计算要补零的字节个数的公式
$\left(3?w+k\right)\mathrm{mod}4=0$
w：是图像的宽度
k：要补几个字节零
mod：取余

#include<stdio.h>
#include<stdlib.h>
#include<Windows.h>
#include<math.h>
#include<string.h>

// 获取文件的后缀名
char* GetFlieExta(char* filename)
{
    int fileLen = strlen(filename);
    int exLen = 0;
    char fileExta[255];
    memset(fileExta, 0, sizeof(fileExta));

    for (int i = fileLen-1; i > 0; i--)
        if (filename[i] == '.'){
            exLen = fileLen - i;
            break;
        }
    strncpy(fileExta, filename + fileLen - exLen, exLen);
    return fileExta;
}


// BGRA颜色结构体
typedef struct tagBGRA
{
	unsigned char blue;          // 该颜色的蓝色分量  (值范围为0-255)
	unsigned char green;         // 该颜色的绿色分量  (值范围为0-255)
	unsigned char red;           // 该颜色的红色分量  (值范围为0-255)
	unsigned char transparency;  // 透明度，在bmp中是保留值，无实际效果
}BGRA, * PBGRA;

// 图像结构体
typedef struct tagIMAGE
{
	unsigned int w;
	unsigned int h;
    BGRA* color;
}IMAGE, * PIMAGE;


// BMP文件的处理

// BMP文件头结构体
typedef struct tagBITMAP_HEAD_INFO
{
    /* bmp文件头的信息，有#的是重点!!*/
    // bmp文件头
    unsigned short  bfType;             // 0x424D，即BM字符串，表明是bmp格式文件
    unsigned int    bfSize;             // ###总的bmp文件大小 以字节为单位     
    unsigned short  bfReserved1;        // 保留，必须设置为0                     
    unsigned short  bfReserved2;        // 保留，必须设置为0 
    unsigned int    bfOffBits;          // ###总的bmp头部的大小（包括位图信息头），即到像素数据的偏移  
    // 位图信息头
    unsigned int    biSize;             // 位图信息头的大小
    unsigned int    biWidth;            // ###图像的宽  
    unsigned int    biHeight;           // ###图像的高  
    unsigned short  biPlanes;           // 颜色平面数，即调色盘数，恒等于1 
    unsigned short  biBitCount;         // ###图片颜色的位数,一般为32
    unsigned int    biCompression;      // 说明图象数据压缩的类型，0为不压缩
    unsigned int    biSizeImage;        // 像素数据所占大小,因为使用BI_RGB，所以设置为0
    unsigned int    biXPelsPerMeter;    // 说明水平分辨率，缺省为0
    unsigned int    biYPelsPerMeter;    // 说明垂直分辨率，缺省为0
    unsigned int    biClrUsed;          // 说明本位图实际使用调色盘的颜色索引数，0表示全部
    unsigned int    biClrImportant;     // 说明本位图重要调色盘的颜色索引数，0表示全都重要
}BITMAP_HEAD_INFO,*PBITMAP_HEAD_INFO;

// 加载BMP图片
IMAGE Image_bmp_load(char* filename)
{
    IMAGE im;
    BITMAP_HEAD_INFO bmpHeadInfo;
    FILE* fp;

    if ((fp = fopen(filename, "rb")) == NULL)
        printf("打开%s文件失败！\n", filename);
    
    // 读取bmp头部
    // bmp文件头
    fread(&bmpHeadInfo.bfType, 1, sizeof(bmpHeadInfo.bfType), fp);
    fread(&bmpHeadInfo.bfSize, 1, sizeof(bmpHeadInfo.bfSize), fp);
    fread(&bmpHeadInfo.bfReserved1, 1, sizeof(bmpHeadInfo.bfReserved1), fp);
    fread(&bmpHeadInfo.bfReserved2, 1, sizeof(bmpHeadInfo.bfReserved2), fp);
    fread(&bmpHeadInfo.bfOffBits, 1, sizeof(bmpHeadInfo.bfOffBits), fp);
    // 位图信息头
    fread(&bmpHeadInfo.biSize, 1, sizeof(bmpHeadInfo.biSize), fp);
    fread(&bmpHeadInfo.biWidth, 1, sizeof(bmpHeadInfo.biWidth), fp);
    fread(&bmpHeadInfo.biHeight, 1, sizeof(bmpHeadInfo.biHeight), fp);
    fread(&bmpHeadInfo.biPlanes, 1, sizeof(bmpHeadInfo.biPlanes), fp);
    fread(&bmpHeadInfo.biBitCount, 1, sizeof(bmpHeadInfo.biBitCount), fp);
    fread(&bmpHeadInfo.biCompression, 1, sizeof(bmpHeadInfo.biCompression), fp);
    fread(&bmpHeadInfo.biSizeImage, 1, sizeof(bmpHeadInfo.biSizeImage), fp);
    fread(&bmpHeadInfo.biXPelsPerMeter, 1, sizeof(bmpHeadInfo.biXPelsPerMeter), fp);
    fread(&bmpHeadInfo.biYPelsPerMeter, 1, sizeof(bmpHeadInfo.biYPelsPerMeter), fp);
    fread(&bmpHeadInfo.biClrUsed, 1, sizeof(bmpHeadInfo.biClrUsed), fp);
    fread(&bmpHeadInfo.biClrImportant, 1, sizeof(bmpHeadInfo.biClrImportant), fp);


    // 读取bmp位图数据
    BGRA* bgra = (BGRA*)malloc(sizeof(BGRA) * (bmpHeadInfo.biWidth * bmpHeadInfo.biHeight));
    fseek(fp, bmpHeadInfo.bfOffBits, SEEK_SET);

    if (bmpHeadInfo.biBitCount == 32)
    {
       for (unsigned int i = 0; i < bmpHeadInfo.biWidth * bmpHeadInfo.biHeight; i++)
            fread(&bgra[i], 1, sizeof(BGRA), fp);
    }
    else if (bmpHeadInfo.biBitCount == 24)
    {
		// 计算每行补几个字节零
         int k = 4 * (3 * bmpHeadInfo.biWidth / 4 + 1) - 3 * bmpHeadInfo.biWidth;
        for (unsigned int i = 0; i < bmpHeadInfo.biWidth * bmpHeadInfo.biHeight; i++)
        {
            if (k != 4 && (ftell(fp)- 54 + k )% (3 * bmpHeadInfo.biWidth + k)==0)
                fseek(fp, ftell(fp) + k, SEEK_SET);

            fread(&bgra[i].blue, 1, sizeof(unsigned char), fp);
            fread(&bgra[i].green, 1, sizeof(unsigned char), fp);
            fread(&bgra[i].red, 1, sizeof(unsigned char), fp);
            bgra[i].transparency = (unsigned char)0xFF;
        }
    }

    im.color = bgra;
    im.w = bmpHeadInfo.biWidth;
    im.h = bmpHeadInfo.biHeight;

    fclose(fp);
    return im;
}

// 保存BMP图片
void Image_bmp_save(char* filename,IMAGE im)
{
    FILE* fp = fopen(filename, "wb");

    unsigned short  bfType = 0x4D42;                // 0x424D，即BM字符串，表明是bmp格式文件
    unsigned int    bfSize = im.w * im.h * 4 + 54;  // ###总的bmp文件大小 以字节为单位     
    unsigned short  bfReserved1 = 0;                // 保留，必须设置为0                     
    unsigned short  bfReserved2 = 0;                // 保留，必须设置为0 
    unsigned int    bfOffBits = 54;                 // ###总的bmp头部的大小（包括位图信息头），即到像素数据的偏移  
    unsigned int    biSize = 40;                    // 位图信息头的大小
    unsigned int    biWidth = im.w;                 // ###图像的宽  
    unsigned int    biHeight = im.h;                // ###图像的高  
    unsigned short  biPlanes = 1;                   // 颜色平面数，即调色盘数，恒等于1 
    unsigned short  biBitCount = 32;                // ###图片颜色的位数,一般为32
    unsigned int    biCompression = 0;              // 说明图象数据压缩的类型，0为不压缩
    unsigned int    biSizeImage = 0;                // 像素数据所占大小,因为使用BI_RGB，所以设置为0
    unsigned int    biXPelsPerMeter = 0;            // 说明水平分辨率，缺省为0
    unsigned int    biYPelsPerMeter = 0;            // 说明垂直分辨率，缺省为0
    unsigned int    biClrUsed = 0;                  // 说明本位图实际使用调色盘的颜色索引数，0表示全部
    unsigned int    biClrImportant = 0;             // 说明本位图重要调色盘的颜色索引数，0表示全都重要

    fwrite(&bfType, 2, 1, fp);
    fwrite(&bfSize, 4, 1, fp);
    fwrite(&bfReserved1, 2, 1, fp);
    fwrite(&bfReserved2, 2, 1, fp);
    fwrite(&bfOffBits, 4, 1, fp);
    fwrite(&biSize, 4, 1, fp);
    fwrite(&biWidth, 4, 1, fp);
    fwrite(&biHeight, 4, 1, fp);
    fwrite(&biPlanes, 2, 1, fp);
    fwrite(&biBitCount, 2, 1, fp);
    fwrite(&biCompression, 4, 1, fp);
    fwrite(&biSizeImage, 4, 1, fp);
    fwrite(&biXPelsPerMeter, 4, 1, fp);
    fwrite(&biYPelsPerMeter, 4, 1, fp);
    fwrite(&biClrUsed, 4, 1, fp);
    fwrite(&biClrImportant, 4, 1, fp);

    fwrite(im.color, sizeof(BGRA) * im.w * im.h, 1, fp);

    fclose(fp);
}


// JPG结构体
// typedef struct tagJPG
// 加载JPG图片
// IMAGE Image_jpg_load(char* filename)
// 保存JPG图片
// void Image_jpg_save(char* filename, IMAGE im)

// 加载图片
IMAGE Image_load(char* filename)
{
    IMAGE im;
    char* fileEx= GetFlieExta(filename);

    if (strcmp(fileEx, ".bmp") == 0)
        im = Image_bmp_load(filename);
    // else if (strcmp(fileEx, ".jpg") == 0)
        // im = Image_jpg_load(filename);

    return im;
}

// 保存图片
void Image_save(char* filename, IMAGE im)
{
    char* fileEx = GetFlieExta(filename);

    if (strcmp(fileEx, ".bmp") == 0)
        Image_bmp_save(filename, im);
    // else if (strcmp(fileEx, ".jpg") == 0)
        // im = Image_jpg_save(filename);
}

// 查看图片
void Image_show(char* filename)
{
    system(filename);
}

// 释放图像结构体
void Image_free(IMAGE im)
{
    free(im.color);
}

// 代码测试
/*
	IMAGE image1 = Image_load("123.bmp");	
	Image_save("01.bmp", image1);
	Image_free(image1);
	Image_show("01.bmp");

*/

二. 图片的形状变换

2.1 图像一维坐标与二维坐标的转换

注意：一维坐标与二维坐标的转换是这篇文章最最基础的核心知识点，贯穿全文的所有代码，一定要理解。图片的核心数据是位图信息头中的RGB值，它可以看作一个巨大的一维数组。图片的第一个像素位于左下角，从左到右，从上到下依次排列。但在图像处理中二维坐标更易处理。所以可将一维数组的下标看作一维坐标，并在需要时转化成二维坐标。

一维坐标转二维坐标
$x=\frac{i}{w}$
$y=i\mathrm{m}\mathrm{o}\mathrm{d}w$
二维坐标转一维坐标
$i=x+y?w$
除法是整除，mod是取余，$i$是一维数组下标，$w$是图片的宽度

2.2 图像的任意大小的缩放

2.2.1 最近邻插值法(不推荐使用)

算法核心：1. 通过新图的坐标计算原图的坐标,逆向映射。
2. 通过原图的坐标计算与该点最邻近的点坐标，原图坐标（0.3,0.3）与它最
邻近的点的坐标为（0,0），即新图坐标中该点的RGB等于原图（0,0）点的RGB。即该点RGB等于最近点的RGB

通过新图的坐标计算原图的坐标
$oldX=newX?\frac{oldW}{newW}$ ???? $oldY=newY?\frac{oldH}{newH}$

// 缩放图片(最近邻插值法)(部分区域有明显的锯齿状，不推荐使用)
IMAGE Transform_shape_nearest(IMAGE im, unsigned int newWidth, unsigned int newHeight)
{
    // 算法核心：通过新图的坐标计算原图的坐标

    float fx, fy;    // f->former  原来（原图坐标）
    float dx, dy;    // d->Decimal 小数（原图坐标的小数部分）

    BGRA* bgra = (BGRA*)malloc(sizeof(BGRA) * newWidth * newHeight);

    // 遍历整张新图
    for (unsigned int i = 0; i < newWidth * newHeight; i++)
    {
        // 通过新图的坐标计算原图的坐标
        fx = (i % newWidth) * ((float)im.w / newWidth);
        fy = (i / newWidth) * ((float)im.h / newHeight);

        // 计算原图坐标的小数部分
        dx = fx - (int)fx;
        dy = fy - (int)fy;

        fx = (dx <= 0.5 ? (int)fx : (int)fx + 1);
        fy = (dy <= 0.5 ? (int)fy : (int)fy + 1);
    
        unsigned int k = fx + fy * im.w; // 注意乘的是原图的w

        if (k >= im.w * im.h)  // 判断k是否越界
            k = im.w * im.h - 1;

        bgra[i].blue = im.color[k].blue;
        bgra[i].green = im.color[k].green;
        bgra[i].red = im.color[k].red;
        bgra[i].transparency = 255;
    }


    free(im.color);
    im.color = bgra;
    im.w = newWidth;
    im.h = newHeight;
    return im;
}
// 代码测试
/*
	IMAGE image1 = Image_load("123.bmp");	
	
	image1 = Transform_shape_nearest(image1, 1200, 600)
	
	Image_free(image1);
	Image_show("01.bmp");
*/

2.2.1 双线性插值法(推荐使用)

算法核心：1. 通过新图的坐标计算原图的坐标，,逆向映射。
2. 通过坐标值，按权重分配像素颜色值,离谁近谁权重大。

计算（0.4,0）
(0.4,0)blue = (0,0)blue $?$ (1-0.4) + (1,0)blue $?$ 0.4
计算（0.4,1）
(0.4,1)blue = (0,1)blue $?$ (1-0.4) + (1,1)blue $?$ 0.4
计算（0.4,0.6）
(0.4,0.6)blue = (0.4,0)blue $?$ (1-0.6) + (0.4,1)blue $?$ 0.6

令四个角上的点的一维坐标为k1,k2,k3,k4。所求点的小数部分为dx,dy可得
$$color=\left(k_1?\left(1?dx\right)+k_2?dx\right)?\left(1?dy\right)+\left(k_3?\left(1?dx\right)+k_4?dx\right)?dy$$

// 缩放图片(双线性插值法)(推荐使用)
IMAGE Transform_shape_linear(IMAGE im, unsigned int newWidth, unsigned int newHeight)
{
    float fx, fy, dx, dy;
    int k1, k2, k3, k4;

    BGRA* bgra = (BGRA*)malloc(sizeof(BGRA) * newWidth * newHeight);

    for (unsigned int i = 0; i < newWidth * newHeight; i++)
    {
        // 通过新图的坐标计算原图的坐标
        fx = (i % newWidth) * ((float)im.w / newWidth);
        fy = (i / newWidth) * ((float)im.h / newHeight);

        dx = fx - (int)fx;
        dy = fy - (int)fy;

        fx = (int)fx;
        fy = (int)fy;

        // 分别计算四个角上点的坐标
        k1 = fx + fy * im.w;
        k2 = fx + 1 + fy * im.w;
        k3 = fx + (fy + 1) * im.w;
        k4 = fx + 1 + (fy + 1) * im.w;

        // 判断是否越界
        if (k1 >= im.w * im.h)
            k1 = im.w * im.h - 1;
        if (k2 >= im.w * im.h)
            k2 = im.w * im.h - 1;
        if (k3 >= im.w * im.h)
            k3 = im.w * im.h - 1;
        if (k4 >= im.w * im.h)
            k4 = im.w * im.h - 1;

        bgra[i].blue = (im.color[k1].blue * (1 - dx) + im.color[k2].blue * dx) * (1 - dy) + (im.color[k3].blue * (1 - dx) + im.color[k4].blue * dx) * dy;
        bgra[i].green = (im.color[k1].green * (1 - dx) + im.color[k2].green * dx) * (1 - dy) + (im.color[k3].green * (1 - dx) + im.color[k4].green * dx) * dy;
        bgra[i].red = (im.color[k1].red * (1 - dx) + im.color[k2].red * dx) * (1 - dy) + (im.color[k3].red * (1 - dx) + im.color[k4].red * dx) * dy;
        bgra[i].transparency = 255;
    }

    free(im.color);
    im.color = bgra;
    im.w = (int)newWidth;
    im.h = (int)newHeight;
    return im;
}

左：最近邻（有明显锯齿），右：双线性
?

2.3 图像的任意角度的旋转

1. 由矩阵旋转可得旋转后的坐标与旋转前坐标的转换

$x^{\prime \prime }=x^{\prime }\cos \theta +y^{\prime }\sin \theta$ ?????$y^{\prime \prime }=y^{\prime }\cos \theta ?x^{\prime }\sin \theta$

2. 由图可知，旋转后的图片离第一象限，$x$轴偏移了$bx$，$y$轴偏移了$by$，所以公式调整为

$x^{\prime \prime }=x^{\prime }\cos \theta +y^{\prime }\sin \theta +bx$?????$y^{\prime \prime }=y^{\prime }\cos \theta ?x^{\prime }\sin \theta +by$

3. 实现的核心：通过原图坐标计算新图坐标，向前映射。
4. 注意旋转后的图片大小也会改变
5. 旋转后的图片会有许多空白区域，最后要对这些空白区域进行填充

// 图像的任意角度的旋转
IMAGE Transform_shape_whirl(IMAGE im, float angle)
{
    // 转角度换成弧度
    angle = 3.141592 * angle / 180;

    float cosnum = (float)cos(angle);   
    float sinnum = (float)sin(angle);

    // 计算原图的四个角的坐标
    int fx1 = 0;
    int fy1 = 0;
    int fx2 = im.w - 1;
    int fy2 = 0;
    int fx3 = 0;
    int fy3 = im.h - 1;
    int fx4 = im.w - 1;
    int fy4 = im.h - 1;

    // 计算旋转后的图像四个角的坐标
    int nx1 = 0;
    int ny1 = 0;
    int nx2 = (int)(fx2 * cosnum + fy2 * sinnum);
    int ny2 = (int)(fy2 * cosnum - fx2 * sinnum);
    int nx3 = (int)(fx3 * cosnum + fy3 * sinnum);
    int ny3 = (int)(fy3 * cosnum - fx3 * sinnum);
    int nx4 = (int)(fx4 * cosnum + fy4 * sinnum);
    int ny4 = (int)(fy4 * cosnum - fx4 * sinnum);

    // 计算旋转后的图像的宽和高
    unsigned int width = abs(max(max(nx1, nx2), max(nx3, nx4))) + abs(min(min(nx1, nx2), min(nx3, nx4))) + 1;
    unsigned int hight = abs(max(max(ny1, ny2), max(ny3, ny4))) + abs(min(min(ny1, ny2), min(ny3, ny4))) + 1;

    // 计算旋转后的图像到第一象限的位置偏移
    int bx = abs(min(min(nx1, nx2), min(nx3, nx4)));
    int by = abs(min(min(ny2, ny3), ny4));

    // 申请并初始化内存空间
    BGRA* bgra = (BGRA*)calloc(width * hight, sizeof(BGRA));
    
    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        // 注意(int)放的位置，不能随便移动,因为y可能为负的小数 
        int k = ((i % im.w) * cosnum + (i / im.w) * sinnum + bx) + (int)((i / im.w) * cosnum - (i % im.w) * sinnum + by) * width;

        // 判断是否越界
        if (k >= width * hight)
            k = width * hight - 1;

        bgra[k].blue = im.color[i].blue;
        bgra[k].green = im.color[i].green;
        bgra[k].red = im.color[i].red;
        bgra[k].transparency = 255;
    }

    // 用邻近的像素填充空白区域
    for (unsigned int i = 0; i < width * hight; i++)
    {
        if (bgra[i].transparency != 255 && bgra[i + 1].transparency == 255)
        {
            bgra[i].blue = bgra[i - 1].blue;
            bgra[i].green = bgra[i - 1].green;
            bgra[i].red = bgra[i - 1].red;
            bgra[i].transparency = 255;
        }
    }

    free(im.color);
    im.color = bgra;
    im.w = width;
    im.h = hight;
    return im;
}

2.4 图像的镜像翻转

图像的镜像翻转分为水平翻转和垂直翻转
水平翻转是$y$坐标不变，$x$坐标翻转
垂直翻转是$x$坐标不变，$y$坐标翻转

#define UPTURN_MODE_HORIZONTAL 0    // 水平翻转
#define UPTURN_MODE_VERTICAL 1      // 垂直翻转
// 图像的镜像翻转
IMAGE Transform_shape_upturn(IMAGE im, int upturn_mode)
{
    BGRA* bgra = (BGRA*)malloc(sizeof(BGRA) * im.w * im.h);

    if(upturn_mode == UPTURN_MODE_HORIZONTAL)
        for (unsigned int i = 0; i < im.w * im.h; i++)
            bgra[i] = im.color[(im.w - 1 - (i % im.w)) + i / im.w * im.w];  // 水平翻转是y坐标不变，x坐标翻转
    else if(upturn_mode == UPTURN_MODE_VERTICAL)
        for (unsigned int i = 0; i < im.w * im.h; i++)
            bgra[i] = im.color[(i % im.w) + (im.h - 1 - i / im.w ) * im.w]; // 垂直翻转是x坐标不变，y坐标翻转

    free(im.color);
    im.color = bgra;
    return im;
}

三. 图片的颜色预处理

图片的颜色预处理中部分算法比较简单所以不做过多的讲解，直接看公式即可明白。

3.1 灰度图（4种算法）

3.1.1 加权法(推荐使用)

$$color=\frac{blue?144+green?587+red?299}{1000}$$

3.1.2 最值法

$$color=max(blue,green,red)$$

3.1.2 均值法

$$color=\frac{blue+green+red}{3}$$

3.1.2 分量法

$$color=(red/green/blue)$$

#define GRAY_MODE_WEIGHT 1           // 加权法(推荐使用)
#define GRAY_MODE_BEST 2             // 最值法
#define GRAY_MODE_AVERAGE 3          // 均值法
#define GRAY_MODE_PART_RED 4         // 分量法_RED
#define GRAY_MODE_PART_GREEN 5       // 分量法_GREEN
#define GRAY_MODE_PART_BLUE 6        // 分量法_BLUE

// 彩色图转灰度图
void Transform_color_grayscale(IMAGE im, int grayscale_mode)
{
    int color = 0;
    switch (grayscale_mode)
    {
    case GRAY_MODE_WEIGHT:
    {
        for (unsigned int i = 0; i < im.w * im.h; i++)
        {
            color = (im.color[i].blue * 114 + im.color[i].green * 587 + im.color[i].red * 299) / 1000;
            im.color[i].blue = color;
            im.color[i].green = color;
            im.color[i].red = color;
        }
        break;
    }

    case GRAY_MODE_BEST:
    {
        for (unsigned int i = 0; i < im.w * im.h; i++)
        {
            color = im.color[i].blue;
            if (color < im.color[i].red)
                color = im.color[i].red;
            if (color < im.color[i].green)
                color = im.color[i].green;
            im.color[i].blue = color;
            im.color[i].green = color;
            im.color[i].red = color;
        }
        break;
    }

    case GRAY_MODE_AVERAGE:
    {
        for (unsigned int i = 0; i < im.w * im.h; i++)
        {
            color = (im.color[i].blue + im.color[i].green + im.color[i].red) / 3;
            im.color[i].blue = color;
            im.color[i].green = color;
            im.color[i].red = color;
        }
        break;
    }

    case GRAY_MODE_PART_RED:
    {
        for (unsigned int i = 0; i < im.w * im.h; i++)
        {
            im.color[i].blue = im.color[i].red;
            im.color[i].green = im.color[i].red;
        }
        break;
    }

    case GRAY_MODE_PART_GREEN:
    {
        for (unsigned int i = 0; i < im.w * im.h; i++)
        {
            im.color[i].blue = im.color[i].green;
            im.color[i].red = im.color[i].green;
        }
        break;
    }

    case GRAY_MODE_PART_BLUE:
    {
        for (unsigned int i = 0; i < im.w * im.h; i++)
        {
            im.color[i].green = im.color[i].blue;
            im.color[i].red = im.color[i].blue;
        }
        break;
    }

    }

}

图片顺序：加权法，最值法，均值法，分量法_RED
?

3.2 二值图

二值图是指由只由黑白两种颜色的像素组成的图片，二值图的质量会直接影响到提取到的图片特征的是否有效，而决定二值图的质量的核心就是阈值的计算，像素值大于等于阈值则取255，小于阈值则取0。

3.2.1 自定义阈值法

// 二值图(自定义阈值法)
void Transform_color_BW_DIY(IMAGE im, unsigned char Threshold)
{
    unsigned char color = 0;
    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        // 先转换成灰度图
        color = (im.color[i].blue * 114 + im.color[i].green * 587 + im.color[i].red * 299) / 1000;
        if (color >= Threshold) // Threshold的值在不同的图片中是不同的
            color = 255;
        else
            color = 0;

        im.color[i].blue = color;
        im.color[i].green = color;
        im.color[i].red = color;
    }
}

3.2.2 最大类间方差法（大津法OTSU）

大津法OTSU是假设阈值将图像分成背景和前景两部分。计算背景和前景之间的类间方差，类间方差越大说明背景和前景之间的差别越大，二值图的效果越好。大津法OSTU，适用双峰直方图。当图像的整体颜色差别不大时，不推荐使用。

类间方差计算公式
$$g=\frac{w_{前}}{1?w_{前}}?{\left(u_{前}?u_{总}\right)}^2$$
$w_{前}$：前景像素数的占比
$u_{前}$：前景的平局灰度
$u_{总}$：图像的平局灰度
$g$：类间方差

整理可得
$$g=\frac{\frac{n_{前}}{n_{总}}}{1?\frac{n_{前}}{n_{总}}}?{\left(\frac{{\sum }_{i=T}^{i=255}{n}_i?i}{n_{总}}?\frac{{\sum }_{i=0}^{i=255}{n}_i?i}{n_{总}}\right)}^2$$
$n_{前}$：前景的像素数
$n_{总}$：总的像素数
$n_i$：对应像素颜色的像素数

这里计算平局灰度，读者可能有点不理解，举个例子：一组7人每人10元，二组6人每人8元，三组11人每人5元，求计算平均每人几元。

解
$$\frac{7?10+6?8+11?5}{7+6+11}=7.2$$

// 二值图(大津法OSTU，适用双峰直方图。当图像的整体颜色差别不大时，不推荐使用)
void Transform_color_BW_OSTU(IMAGE im)
{
    // 公式：g = w0 / (1 - w0) * (u0 - u)* (u0 - u) 当g最大时取到阈值T
    int colorMap[256] = { 0 };
    float w0 = 0; // 前景像素数的占比
    unsigned int u0 = 0; // 前景的平局灰度（灰度值*其对应的素数个数）的累加/前景的像素个数
    unsigned int u = 0;  // 图像的平局灰度（灰度值*其对应的素数个数）的累加/总的像素个数
    float g = 0;  // 方差
    unsigned char T = 0;    // 阈值

    // 创建灰度直方图
    for (unsigned int i = 0; i < im.w * im.h; i++)
        colorMap[im.color[i].blue] += 1;
  
    for (int i = 0; i < 256; i++)
         u += colorMap[i] * i;  // u暂时计算累加
    u /= (im.h * im.w);

    // 遍历 0-255 寻找合适的阈值 
    for(unsigned int m = 0 ; m < 256; m++)
    {
        for (int n = m; n < 256; n++)
        {
            w0 += colorMap[n]; // w0暂时计算，保存前景的所有像素个数
            u0 += colorMap[n] * n; // u0暂时计算累加
        }
        u0 /= w0;
        w0 /= (im.h * im.w);
        
        if((w0 / (1 - w0) * (u0 - u)* (u0 - u)) > g)
            g = w0 / (1 - w0) * (u0 - u)* (u0 - u), T = m;
    }
 Transform_color_BW_DIY(im, T);
}

3.2.3 三角法TRIANGLE

如图，在灰度直方图中，找到最低点和最高点，计算出这两点的函数（黄色直线），分别计算最低点和最高点之间的点到黄色直线的距离$d$，当$d$取到最大时，取到$T$阈值。三角法TRIANGLE，适用单峰直方图。当图像的整体颜色差别不大时，不推荐使用

计算公式
? ? ? ? $d=\frac{|A{x}_0+B{y}_0+C|}{\sqrt{A^2+B^2}}$ ?? ? ? $y=kx+b$
联立可得
$$d=\frac{\left|?k{x}_0+y_0?b\right|}{\sqrt{1+k^2}}$$

// 二值图(三角法TRIANGLE，适用单峰直方图。当图像的整体颜色差别不大时，不推荐使用)
void Transform_color_BW_TRIANGLE(IMAGE im)
{
    int colorMap[256] = { 0 };
    unsigned char minColor = 0;
    unsigned int minCount = 0;
    unsigned char maxColor = 0;
    unsigned int maxCount = 0;
    int d = 0;  // 最短距离
    unsigned char T = 0; // 阈值

    // 创建灰度直方图
    for (unsigned int i = 0; i < im.w * im.h; i++)
        colorMap[im.color[i].blue] += 1;

    for (int i = 0; i < 256; i++)
    {
        if (im.color[i].blue < minColor)
            minColor = im.color[i].blue, minCount = colorMap[im.color[i].blue];
        if (im.color[i].blue > maxColor)
            maxColor = im.color[i].blue, maxCount = colorMap[im.color[i].blue];
    }

    float k = ((float)maxCount - minCount) / ((float)maxColor - minColor);
    float b = maxCount - k * maxColor;
   
    // 遍历寻找最近距离 
    for (unsigned int n = minColor; n <= maxColor; n++)
        if (abs((int)(-k * n + colorMap[n] - b)) / sqrt((double)(1 + k * k)) > b)
            b = abs((int)(-k * n + colorMap[n] - b)) / sqrt((double)(1 + k * k)), T = n;

    Transform_color_BW_DIY(im, T);
}

识别出了部分云彩

3.2.4 自适应阈值法

自适应阈值法是一种，类卷积的操作（卷积会在下文详细讲解），其核心是选取一定大小的区域(区域的边长为奇数)，用这区域的像素平均值来作为这一个区域的阈值。自适应阈值法，当图片线条多且密时，不推荐使用

// 二值图(自适应阈值法，areaSize=25较合适，当图片线条多且密时，不推荐使用)
IMAGE Transform_color_BW_Adaptive(IMAGE im, int areaSize)
{
    // areaSize为区域的大小，区域越大，效果图的细节越好，areaSize=25较合适
    BGRA* bgra = (BGRA*)malloc(sizeof(BGRA) * im.w * im.h);
    int* p = (int*)malloc(sizeof(int) * areaSize); // p->position 位置坐标
    int k = (int)(sqrt((double)areaSize)) / 2;  // 重合区域边长的一半

    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        // 计算与卷积和对应重合区域的坐标
        int t = 0; // 记录p的下标
        for (int n = k; n >= -k; n--)
            for (int m = -k; m <= k; m++)
            {
                p[t] = ((i % im.w) + m) + (i / im.w + n) * im.w;
                t++;
            }

        // 判断是否越界
        for (int j = 0; j < areaSize; j++)
            if (p[j] < 0 || p[j] >= im.w * im.h)
                p[j] = i;

        unsigned int color = 0;
        for (int j = 0; j < areaSize; j++)
            color += im.color[p[j]].blue;
        color /= areaSize;

        if (im.color[i].blue >= color)
            bgra[i].blue = 255;
        else
            bgra[i].blue = 0;

        bgra[i].green = bgra[i].blue;
        bgra[i].red = bgra[i].blue;
    }

    free(p);
    free(im.color);
    im.color = bgra;
    return im;
}

处理山

处理人脸

3.2.5 拓展：用二值图表示灰度变化

核心思想是通过像素的疏密来表示灰度的变化，其实现的代码和自适应阈值法的代码极为相识

// 二值图(用二值图表示灰度变化,areaSize=25较合适)
void Transform_color_BW_grayscale(IMAGE im, int areaSize)
{
	// areaSize为区域的大小，区域越大，效果图的细节越好，areaSize=25较合适
    int* p = (int*)malloc(sizeof(int) * areaSize); // p->position 位置坐标
    int k = (int)(sqrt((double)areaSize)) / 2;  // 重合区域边长的一半

    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        // 计算与卷积和对应重合区域的坐标
        int t = 0; // 记录p的下标
        for (int n = k; n >= -k; n--)
            for (int m = -k; m <= k; m++)
            {
                p[t] = ((i % im.w) + m) + (i / im.w + n) * im.w;
                t++;
            }

        // 判断是否越界
        for (int j = 0; j < areaSize; j++)
            if (p[j] < 0 || p[j] >= im.w * im.h)
                p[j] = i;

        unsigned int color = 0;
        for (int j = 0; j < areaSize; j++)
            color += im.color[p[j]].blue;
        color /= areaSize;

        if (im.color[i].blue >= color)
            im.color[i].blue = 255;
        else
            im.color[i].blue = 0;

        im.color[i].green = im.color[i].blue;
        im.color[i].red = im.color[i].blue;
    }

    free(p);
}

这是用彩色图转的看起来效果较好，也可以用灰度图

3.3 反色

$colo{r}^{\prime \prime }=255?colo{r}^{\prime }$

// 反色
void Transform_color_opposite(IMAGE im)    
{
    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        im.color[i].green = 255 - im.color[i].green;
        im.color[i].blue = 255 - im.color[i].blue;
        im.color[i].red = 255 - im.color[i].red;
    }

}

3.4 直方图均衡化

1. 直方图均衡化：是一种增强图像对比度的方法，即让亮的地方更亮，暗的地方更暗，主要用于图像整体颜色偏淡，图像整体太亮，图像整体太暗，等图像整体变化度不够的情况。由图可知，其主要思想是将一副图像的颜色直方图分布变成近似均匀分布，从而增强图像的对比度。
2. 直方图均衡化核心公式

$$color=\frac{\left(RGB\max ?RGB\min \right)}{N}\sum _{i=0}^{k}ni$$
$RGB\max$是：最大的RGB值一般为255
$RGB\min$是：最小的RGB值一般为0
N是：总的像素数一般为$w?h$
$k$ 是：color的RGB值
${\sum }_{i=0}^{k}ni$ 是：所有小于color的RGB值的像素的个数之和

整理后的公式

$$color=\frac{255?0}{w?h}\sum _{i=0}^{k}ni$$

3.两种不同的直方图均衡化（实现代码差不多，但效果区别很大）

// 直方图均衡化（分步计算，效果更加柔和）
void Transform_color_Histogram_part(IMAGE im)
{
    // 公式：均衡后的颜色值=（最大颜色位255-最小颜色位0）*小于等于该颜色值的像素数量的累加/图片总的像素数

    int Accumulate = 0;             // 保存累加的值
    unsigned char color = 0;                  // 保存颜色的值
    int allBlue[256] = { 0 };       // 保存蓝色直方图
    int allGreen[256] = { 0 };      // 保存绿色直方图
    int allRed[256] = { 0 };        // 保存红色直方图

    // 数组下标等于RBG值极大的简化了计算
    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        allBlue[im.color[i].blue] += 1;
        allGreen[im.color[i].green] += 1;
        allRed[im.color[i].red] += 1;
    }

    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        //blue
        for (int j = 0; j <= im.color[i].blue; j++) // 累加计算
            Accumulate += allBlue[j];
        color = (255 - 0) * Accumulate / (im.w * im.h);
        im.color[i].blue = color;
        Accumulate = 0;

        //green
        for (int j = 0; j <= im.color[i].green; j++)// 累加计算
            Accumulate += allGreen[j];
        color = (255 - 0) * Accumulate / (im.w * im.h);
        im.color[i].green = color;
        Accumulate = 0;

        //red
        for (int j = 0; j <= im.color[i].red; j++)// 累加计算
            Accumulate += allRed[j];
        color = (255 - 0) * Accumulate / (im.w * im.h);
        im.color[i].red = color;
        Accumulate = 0;
    }
}

// 直方图均衡化（整体计算，效果更加尖锐）
void Transform_color_Histogram_all(IMAGE im)
{    
    int Accumulate = 0;         // 保存累加的值
    unsigned char color = 0;              // 保存颜色的值
    int allColor[256] = { 0 };  // 保存所有颜色直方图

    // 数组下标等于RBG值极大的简化了计算
    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        allColor[im.color[i].blue] += 1;
        allColor[im.color[i].green] += 1;
        allColor[im.color[i].red] += 1;
    }

    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        //blue
        for (int j = 0; j <= im.color[i].blue; j++)// 累加计算
            Accumulate += allColor[j];
        color = (255 - 0) * Accumulate / (im.w * im.h);
        im.color[i].blue = color;
        Accumulate = 0;

        //green
        for (int j = 0; j <= im.color[i].green; j++)// 累加计算
            Accumulate += allColor[j];
        color = (255 - 0) * Accumulate / (im.w * im.h);
        im.color[i].green = color;
        Accumulate = 0;

        //red
        for (int j = 0; j <= im.color[i].red; j++)// 累加计算
            Accumulate += allColor[j];
        color = (255 - 0) * Accumulate / (im.w * im.h);
        im.color[i].red = color;
        Accumulate = 0;
    }
}

分步计算，效果更加柔和

整体计算，效果更加尖锐

注意：上文写二值图时介绍的大津法和三角法，都在图像的整体颜色差别不大时处理效果不佳，而直方均衡化就抑制这问题。

四. 图像的卷积操作

4.1 卷积核

1. 卷积核是图像处理时，给定输入图像，输入图像中一个小区域中像素加权平均后成为输出图像中的每个对应像素，其中权值由一个函数定义，这个函数称为卷积核。通俗的讲就是一个一维或二维数组，但注意卷积核一般是$3?3$，$5?5$，$7?7$ 奇数大小的。卷积核中的数加起来一般为0。 灵活运用不同的卷积核可以实现很多功能。
   
2. 卷积核的操作和计算方法，如图所示，$3?3$的卷积核在$w?h$的$6?6$图像上进行移动，将卷积核与图像重合部分的数，分别相乘再相加，即为卷积核中心的值（黄色框区域）当该值小于0时，一般会取绝对值。
3. 不与卷积核重合的部分，即超出图像边界的部分的值一般为设为0，或等于卷积核中心的值。
4. 类卷积操作，指并没有用到卷积核，但操作方法与卷积操作类似

4.1.1 自定义卷积操作

卷积操作的代码实现具有很高的重复性，所以先写一个自定义卷积操作（代码与上文二值图的自定义阈值法的实现相类似）可以极大的缩减代码。类卷积操作的代码实现还是要手写。

// 判断像素值的范围
unsigned char Tool_RBG(int BRRA)
{
    if (BRRA > 255)
        return (unsigned char)255;
    else if (BRRA < 0)
        return (unsigned char)0;
    else
        return (unsigned char)BRRA;
}

// 卷积操作（自定义）
IMAGE Kernels_use_DIY(IMAGE im, double* kernels, int areaSize, int modulus)
{
	// kernels卷积核
	// areaSize区域的大小
	// modulus最后乘的系数
	
    BGRA* bgra = (BGRA*)malloc(sizeof(BGRA) * im.w * im.h);
    int* p = (int*)malloc(sizeof(int) * areaSize); // p->position 位置坐标
    int k = (int)(sqrt((double)areaSize)) / 2;  // 重合区域边长的一半

    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        // 计算与卷积和对应重合区域的坐标
        int t = 0; // 记录p的下标
        for(int n = k; n >= -k; n--)
            for (int m = -k; m <= k; m++)
                p[t] = ((i % im.w) + m) + (i / im.w + n) * im.w, t++;
                
        // 判断是否越界
        for (int j = 0; j < areaSize; j++) 
            if (p[j] < 0 || p[j] >= im.w * im.h)
                p[j] = i;

		// 相乘相加
        int blue = 0, green = 0, red = 0;
        for (int j = 0; j < areaSize; j++)
        {
            blue += im.color[p[j]].blue * kernels[j];
            green += im.color[p[j]].green * kernels[j];
            red += im.color[p[j]].red * kernels[j];
        }
        
        bgra[i].blue = Tool_RBG(blue * modulus);
        bgra[i].green = Tool_RBG(green * modulus);
        bgra[i].red = Tool_RBG(red * modulus);
    }

    free(p);
    free(im.color);
    im.color = bgra;
    return im;
}

4.2 滤波

图像滤波即在尽量保留图像细节特征的条件下对目标图像的噪声进行抑制。

4.2.1 中值滤波

中值滤波主要用于处理脉冲噪声（椒盐噪声），它是属于随机噪声，即在图像上随机出现的黑点或白点。中值滤波，是将卷积区域内的值按顺序排列，并取中值，作为该点的像素。

//中值滤波
IMAGE Wavefiltering_Median(IMAGE im)
{
    BGRA* bgra = (BGRA*)malloc(sizeof(BGRA) * im.w * im.h);

    for (int unsigned i = 0; i < im.w * im.h; i++)
    {
        // 与卷积和对应重合区域的坐标
        int p[9] =   // p->position 位置坐标
        {
            i + im.w - 1,i + im.w,i + im.w + 1,
            i - 1,i,i + 1,
            i - im.w - 1,i - im.w,i - im.w + 1
        };

        // 判断是否越界
        for (int j = 0; j < 9; j++)
            if (p[j] < 0 || p[j] >= im.w * im.h)
                p[j] = i;

        // 取颜色
        int color[9] =
        {
            im.color[p[0]].blue, im.color[p[1]].blue, im.color[p[2]].blue,
            im.color[p[3]].blue, im.color[p[4]].blue, im.color[p[5]].blue,
            im.color[p[6]].blue, im.color[p[7]].blue, im.color[p[8]].blue
        };

        for (int n = 0; n < 9; n++)  // 对颜色进行排序
            for (int m = n; m < 9; m++)
                if (color[n] > color[m])   // 异或交换不能用等于号
                {
                    color[n] ^= color[m];
                    color[m] ^= color[n];
                    color[n] ^= color[m];
                }

        bgra[i].blue = color[4];
        bgra[i].green = color[4];
        bgra[i].red = color[4];
    }

    free(im.color);
    im.color = bgra;
    return im;
}

滤波前后比较

4.2.2 高斯滤波

高斯噪声是指它的概率密度函数服从高斯分布（即正态分布）的一类噪声。常见的高斯噪声包括起伏噪声、宇宙噪声、热噪声和散粒噪声等等。

高斯滤波公式
$$\frac{1}{16}?\begin{array}{ccc}1& 2& 1\\ & & \\ 2& 4& 2\\ & & \\ 1& 2& 1\\ & & \end{array}?$$

//高斯滤波卷积核
double KERNELS_Wave_Gauss[9] =
{
    1, 2, 1,
    2, 4, 2,
    1, 2 ,1
};

//高斯滤波
IMAGE Wavefiltering_Gauss(IMAGE im)
{
    im = Kernels_use_DIY(im, KERNELS_Wave_Gauss, 9, 1.0 / 16);
    return im;
}

效果和中值滤波差不多

4.2.3 低通滤波

低通滤波是增强图像的低频成分，阻挡高频成分，使图像变得平滑与模糊。

// 低通滤波卷积核 LP1
double KERNELS_Wave_LowPass_LP1[9] =
{
    1 / 9.0, 1 / 9.0, 1 / 9.0,
    1 / 9.0, 1 / 9.0, 1 / 9.0,
    1 / 9.0, 1 / 9.0, 1 / 9.0
};

// 低通滤波卷积核 LP2
double KERNELS_Wave_LowPass_LP2[9] =
{
    1 / 10.0, 1 / 10.0, 1 / 10.0,
    1 / 10.0, 1 / 5.0, 1 / 10.0,
    1 / 10.0, 1 / 10.0, 1 / 10.0
};

// 低通滤波卷积核 LP3
double KERNELS_Wave_LowPass_LP3[9] =
{
    1 / 16.0, 1 / 8.0, 1 / 16.0,
    1 / 8.0, 1 / 4.0, 1 / 8.0,
    1 / 16.0, 1 / 8.0, 1 / 16.0
};


// 低通滤波
IMAGE Wavefiltering_LowPass(IMAGE im, double* kernels)
{
    im = Kernels_use_DIY(im, kernels, 9, 1);
    return im;
}

LP1

LP2

LP3

4.2.4 高通滤波

高通滤波是增强图像的高频成分，阻挡低频成分，使图像变得锐化与清晰。

// 高通滤波卷积核 HP1
double KERNELS_Wave_HighPass_HP1[9] =
{
    -1, -1, -1,
    -1, 9, -1,
    -1, -1 ,-1
};

// 高通滤波卷积核 HP2
double KERNELS_Wave_HighPass_HP2[9] =
{
    0, -1, 0,
   -1, 5, -1,
    0, -1 ,0
};

// 高通滤波卷积核 HP3
double KERNELS_Wave_HighPass_HP3[9] =
{
    1, -2, 1,
   -2, 5, -2,
    1, -2 ,1
};

// 高通滤波
IMAGE Wavefiltering_HighPass(IMAGE im, double* kernels)
{
    im = Kernels_use_DIY(im, kernels, 9, 1);
    return im;
}

HP1

HP2

HP3

4.2.5 均值滤波

和低通滤波差不多，都是使图像变得平滑与模糊。

// 均值滤波卷积核
double KERNELS_Wave_Average[25] =
{
  1, 1, 1, 1, 1,
  1, 1, 1, 1, 1,
  1, 1, 1, 1, 1,
  1, 1, 1, 1, 1,
  1, 1, 1, 1, 1
};

// 均值滤波
IMAGE Wavefiltering_Average(IMAGE im)
{
    im = Kernels_use_DIY(im, KERNELS_Wave_Average, 25, 1.0 / 25);
    return im;
}

4.3 边缘检测

边缘检测就是找出图像中像素亮度发生剧烈变化的区域，通常这些区域表现出就是图像的轮廓。

4.3.1 差分边缘检测

差分边缘检测有三个卷积核，分别用于检测水平方向，垂直方向，垂直和水平方向。

// 差分垂直边缘检测卷积核
double KERNELS_Edge_difference_vertical[9] =
{
    0, 0, 0,
   -1, 1, 0,
    0, 0, 0
};

// 差分水平边缘检测卷积核
double KERNELS_Edge_difference_horizontal[9] =
{
    0,-1, 0,
    0, 1, 0,
    0, 0, 0
};

// 差分垂直和水平边缘检测卷积核
double KERNELS_Edge_difference_VH[9] =
{
   -1, 0, 0,
    0, 1, 0,
    0, 0, 0
};


// 差分边缘检测
IMAGE Edge_detection_difference(IMAGE im, double* kernels)
{
    im = Kernels_use_DIY(im, kernels, 9, 1);
    return im;
}

差分垂直和水平边缘检测

4.3.2 Sobel边缘检测

差分边缘检测有两个卷积核，分别对应XY方向，而该点的像素值为

$$color=\sqrt{x^2+y^2}$$

// Sobel X边缘检测卷积核
double KERNELS_Edge_Sobel_X[9] =
{
    -1, 0, 1,
   - 2, 0, 2,
    -1, 0, 1
};

// Sobel Y边缘检测卷积核
double KERNELS_Edge_Sobel_Y[9] =
{
   -1, -2, -1,
    0, 0, 0,
    1, 2, 1
};

// Sobel边缘检测
IMAGE Kernels_use_Edge_Sobel(IMAGE im, double* kernels1, double* kernels2)
{
    BGRA* bgra = (BGRA*)malloc(sizeof(BGRA) * im.w * im.h);

    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        // 与卷积和对应重合区域的坐标
        int p[9] =   // p->position 位置坐标
        {
            i + im.w - 1,i + im.w,i + im.w + 1,
            i - 1,i,i + 1,
            i - im.w - 1,i - im.w,i - im.w + 1
        };

        for (int j = 0; j < 9; j++) // 判断是否越界
            if (p[j] < 0 || p[j] >= im.w * im.h)
                p[j] = i;


        unsigned char color1 = Tool_RBG(im.color[p[0]].blue * kernels1[0] + im.color[p[1]].blue * kernels1[1] + im.color[p[2]].blue * kernels1[2] + im.color[p[3]].blue * kernels1[3] + im.color[p[4]].blue * kernels1[4] + im.color[p[5]].blue * kernels1[5] + im.color[p[6]].blue * kernels1[6] + im.color[p[7]].blue * kernels1[7] + im.color[p[8]].blue * kernels1[8]);

        if (kernels2 != NULL)
        {
            unsigned char  color2 = Tool_RBG(im.color[p[0]].blue * kernels2[0] + im.color[p[1]].blue * kernels2[1] + im.color[p[2]].blue * kernels2[2] + im.color[p[3]].blue * kernels2[3] + im.color[p[4]].blue * kernels2[4] + im.color[p[5]].blue * kernels2[5] + im.color[p[6]].blue * kernels2[6] + im.color[p[7]].blue * kernels2[7] + im.color[p[8]].blue * kernels2[8]);
            color1 = Tool_RBG(sqrt(color1 * color1 + color2 * color2));
        }

        bgra[i].blue = color1;
        bgra[i].green = color1;
        bgra[i].red = color1;
    }

    free(im.color);
    im.color = bgra;
    return im;
}

X方向

Y方向

XY方向

4.3.3 Laplace边缘检测

Laplace边缘检测的最大特点是它的边缘检测与方向无关。

// Laplace边缘检测卷积核 LAP1
double KERNELS_Edge_Laplace_LAP1[9] =
{
    0, 1, 0,
    1, -4, 1,
    0, 1, 0
};

// Laplace边缘检测卷积核 LAP2
double KERNELS_Edge_Laplace_LAP2[9] =
{
   -1, -1, -1,
   -1, 8, -1,
   -1, -1, -1
};

// Laplace边缘检测卷积核 LAP3
double KERNELS_Edge_Laplace_LAP3[9] =
{
   -1, -1, -1,
   -1, 9, -1,
   -1, -1, -1
};

// Laplace边缘检测卷积核 LAP4
double KERNELS_Edge_Laplace_LAP4[9] =
{
    1, -2, 1,
   -2, 8, -2,
    1, -2, 1
};


// Laplace边缘检测 
IMAGE Edge_detection_Laplace(IMAGE im, double* kernels)
{
    im = Kernels_use_DIY(im, kernels, 9, 1);
    return im;
}

LAP1

LAP2

LAP3

LAP4

4.4 形态学图像处理

1. 形态学图像处理包括，图形腐蚀，膨胀，开闭运算等，一般都是对二值图进行处理。且二值图是黑底白字（可通过反色转换）
2. 形态学图像处理中的用的卷积核种类很多，可根据实际自己设计。

4.4.1 腐蚀

1. 图像腐蚀是将图像中的高亮区域即白色部分进行缩减细化，其运行结果图比原图的高亮区域更小。即图像轮廓变细了。
2. 腐蚀用的卷积核，用的是十字形

// 腐蚀卷积核
double KERNELS_Morphology_Erosion_cross[9] =
{
    0, 1, 0,
    1, 1, 1,
    0, 1, 0
};

3. 腐蚀卷积核的使用规则，除去卷积核的中心，将其它8个数与图像上对应的数相乘再相加，如果结果小于255乘卷积核中8个数的和的积，那么就将卷积核的中心所对应的图像上的数变为0 。(即边上四个1中有一个所对应的数是0就成立)
4. 腐蚀用的公式

$$\sum _{n=1}^{8}kn?pn<255?\sum _{n=1}^{8}kn$$
kn是：卷积核中的数
pn是：卷积核所对应的图像上的数
注意：这八个数字不包含卷积核的中心

十字形腐蚀卷积核的公式
$$\left(k_1?p_1+k_2?p_2+k_3?p_3+k_4?p_4\right)<255?\left(k_1+k_2+k_3+k_4\right)$$
注意：这里是方便读者理解，但实际代码实现中为了代码的通用性，是用原公式实现的。

// 腐蚀
IMAGE Morphology_Erosion(IMAGE im, double* kernels)
{
    BGRA* bgra = (BGRA*)malloc(sizeof(BGRA) * im.w * im.h);

    for (unsigned int i = 0; i < im.w * im.h; i++)
    {
        // 与卷积和对应重合区域的坐标
        int p[9] =   // p->position 位置坐标
        {
            i + im.w - 1,i + im.w,i + im.w + 1,
            i - 1,i,i + 1,
            i - im.w - 1,i - im.w,i - im.w + 1
        };

        for (int j = 0; j < 9; j++) // 判断是否越界
            if (p[j] < 0 || p[j] >= im.w * im.h)
                p[j] = i;

        // 判断是否腐蚀（式子很长，但简单）
        if ((im.color[i].blue == 255) && (im.color[p[0]].blue * kernels[0] + im.color[p[1]].blue * kernels[1] + im.color[p[2]].blue * kernels[2] + im.color[p[3]].blue * kernels[3] + im.color[p[5]].blue * kernels[5] + im.color[p[6]].blue * kernels[6] + im.color[p[7]].blue * kernels[7] + im.color[p[8]].blue * kernels[8]) < 255 * (kernels[0] + kernels[1] + kernels[2] + kernels[3] + kernels[5] + kernels[6] + kernels[7] + kernels[8]))
        {
            bgra[i].blue = 0;
            bgra[i].green = 0;
            bgra[i].red = 0;
        }
        else
        {
            bgra[i].blue = im.color[i].blue;
            bgra[i].green = im.color[i].green;
            bgra[i].red = im.color[i].red;
        }
    }

    free(im.color);
    im.color = bgra;
    return im;
}

腐蚀前后效果比较

4.4.2 膨胀

1. 图像膨胀是将图像中的高亮区域即白色部分进行扩张添加像素值，其运行结果图比原图的高亮区域更大。即图像轮廓变粗了。
2. 膨胀用的卷积核与腐蚀用的一样，都是十字形。

// 膨胀卷积核
double KERNELS_Morphology_Dilation_cross[9] =
{
    0, 1, 0,
    1, 1, 1,
    0, 1, 0
};