图像编码之赫夫曼编码
图像编码
写在前面,这一次的实验相当有难度,完全自己写体会很深刻,把C/C++对于文件读写部分的内容相当全面的覆盖到了,并且也对算法的设计有较高要求。
按照惯例我们还是先说说为什么要对图像进行编码。现在4K逐渐开始普及了,8K的产品也开始销售了,但是这些超高清背后的数据量,你有考虑过嘛?
我们拿现在已经不算稀奇的FullHD(1080p)来做一笔计算:
一幅(帧)图像的字节数:
120分钟的电影:
如果我们用比较好的50G蓝光光盘保存,需要:
张蓝光光碟。
可见,除非Sony宣布蓝光光盘免费,否则看一场电影的门槛可能都要上百了。正因如此,我们才要将图像进行“压缩”,图像能够被压缩的客观条件,就是“冗余”的存在。
一幅图片从不同角度来看会存在多种多样的冗余信息(可以参考这篇科普),我们适当的删减掉不容易被人眼注意到的部分,就可以在不减少信息量的情况下压缩图片的大小。当然这个压缩的具体效果肯定还是与删去的冗余信息量有关的,一幅JPEG压缩的图像在多数情况下都会比原图片看起来“糊”不少,就是因为JPEG的压缩比能够达到26:1左右。
图像编码有非常多的方式,这次的实验内容,就是完成一种基于熵编码的算法——赫夫曼编码。原理这方面不再赘述了,大把的博客以及书籍可供参考,我们直接进入实验内容。
一、实验内容
将一幅给定的灰度图像进行Huffman编解码
众所周知,字数越少,事情越大,这短短一行的实验要求耗费了我三整天的时间。根据我做这个实验过程中遇到的问题,我来将这个实验内容拆分一下:
- 编写函数,获得赫夫曼编码表
- 根据赫夫曼编码,将像素数据用编码表示
- 将编码后的图片数据写入文件
- 读取压缩后的图像文件和赫夫曼编码表
- 利用赫夫曼编码表还原图像文件
二、代码实现与分析
1、编写函数,获得赫夫曼编码表
别看这是唯一一个和标题名字重合的步骤,实际上却是最简单的,当然这也是因为我上学期学习数据结构时就已经完成了赫夫曼树、赫夫曼编码的函数,这里我只进行了简单的修改,便得到了码表。
简单的讲一下代码,为什么这么写可以参考我老师的mooc
typedef struct {
float weight; //结点权值
int parent, lchild, rchild; //指向双亲和左右孩子的指针
}HTNode,*HuffmanTree; //动态分配数组赫夫曼树
typedef char **HuffmanCode; //动态分配数组存储赫夫曼编码
typedef struct {
HuffmanCode code; //码字数组,从下标1开始
int *length; //码长数组,从下标1开始
}HuffmanCodeTable; //赫夫曼编码表
上面是结构体的定义,我们的思路是先构造赫夫曼树,再根据赫夫曼树从根节点开始走到各个结点的路径得到赫夫曼编码(左0右1)。上这门课的时候要求使用C实现,实际上用C++的string类型可能更加方便。
void HuffmanCoding(HuffmanTree &HT, HuffmanCode &HC, HuffmanCodeTable &HCT, float *w, int n)
{ //w存放n个字符的权值,构造赫夫曼树HT,并求出w中n个数据对应的赫夫曼编码表HCT
int i, j, m, s1, s2;
if (n <= 1)
return;
HTNode *p; //指针
m = 2 * n - 1; //m是总的结点数
HT = (HuffmanTree)malloc((m + 1) * sizeof(HTNode)); //0号单元作为标记
for (p = HT + 1, i = 1; i <= n; i++, p++, w++)
{
*p = { *w,0,0,0 }; //叶子结点初始化
}
for (i = n + 1; i <= m; i++, p++)
{
*p = { 0,0,0,0 }; //非叶子结点初始化
}
for (i = n + 1; i <= m; i++)
{
Select(HT, i - 1, s1, s2); //寻找最小值下标s1和次小值下标s2
HT[s1].parent = i;
HT[s2].parent = i;
HT[i].lchild = s1;
HT[i].rchild = s2;
HT[i].weight = HT[s1].weight + HT[s2].weight;
} //赫夫曼树建立完毕
//求赫夫曼编码的准备工作
HC = (HuffmanCode)malloc((n+1) * sizeof(char*)); //分配n个字符串编码的头指针向量
char *cd = (char*)malloc(n * sizeof(char)); //分配求当前编码的工作空间
cd[n - 1] = '\0'; //从右向左逐位存放编码,首先存放编码结束符
int start, c, f;
for (i = 1; i <= n; i++) //按照w中的顺序给出赫夫曼编码
{
start = n - 1; //初始化编码起始指针
for (c = i, f = HT[i].parent; f != 0; c = f, f = HT[f].parent) //从叶子结点到根结点逆向求编码
{
if (HT[f].lchild == c)
cd[--start] = '0'; //左分支标0
else
cd[--start] = '1'; //右分支标1
}
HC[i] = (char*)malloc((n - start) * sizeof(char)); //为第i个编码分配空间,(n-start)为编码长度
strcpy(HC[i], &cd[start]); //从cd复制编码串到HC
}
free(cd);
//对赫夫曼编码的码字求长度,得到赫夫曼编码表
HCT.code = HC;
HCT.length = (int*)malloc(sizeof(int)*(n+1));
int count;
for (i = 1; i <= n; i++)
{
count = 0; //实际上可以用strlen()
for (j = 0; HC[i][j] != '\0'; j++)
{
count += 1;
}
HCT.length[i] = count;
}
}
赫夫曼编码中有一个Select函数是用来求所有没有父结点的结点中权值最小的两个,从而构成新的结点。我的实现写复杂了,想参考可以见文末传送门,建议写一个排序来实现,效率会更高。
2、根据赫夫曼编码,将像素数据用编码表示
上面的函数只是能够实现赫夫曼编码,实际上我们还没有开始调用,下面来说一下对图像调用赫夫曼编码的想法。我们压缩压的是像素数据域的内容,显然,如果我们把出现次数多的像素值赋予短编码,出现次数少的像素值赋予长编码,那么从均值的角度看,我们还是成功的减小了总长度,也就实现了压缩的目的。
因此,我们应当统计图像中各个像素值(我们这里只考虑灰度图像,故为灰度值)出现的频率,然后对频率进行赫夫曼编码。
for (i = 0; i < src.imagew; i++) //统计像素出现次数
for (j = 0; j < src.imageh; j++)
{
count_pixel[src.pDataAt(i)[j]] += 1;
}
for (i = 0; i < 256; i++) //计算像素出现频率
{
p[i] = (float)count_pixel[i] / (float)n;
}
HuffmanCoding(HT, HC, HCT, p, 256); //赫夫曼编码函数
这样我们就成功的获得了HCT,赫夫曼编码表,它由一个码字数组和一个码长数组构成(均从下标1开始),将其写入文件后,如下:(使用'\t'
分隔码字与码长)
下面我们用它来将像素数据域的内容编码。不过在转换编码前,我们应该先做一个思考,应当用什么格式来存储编码后的内容呢?
首先我们要明确,写入文件一定是以二进制形式完成的,因为ASCII码表示的'0'
和'1'
本身就占一个字节,如果以文本格式输出则完全是与压缩的初衷背向而驰。所以,我们要选一个很容易转换为二进制的格式,字符数组便成了首选。
其次,char*
也可以表示字符数组,string
也可以表示字符数组,具体用哪个更好呢。我们需要将所有像素编码后的结果合并到一个大字符串中,因此对于char*
而言,我们要用strcat(des,src)
,对于string
而言,我们可以用s1+=s2
。经过我的实测,strcat执行时间超过40s,string的+=只需要2s,这么大的差异也确实出乎我的意料,不过这篇博客中的测试也显示了strcat的效率不尽人意。
string data;
for (i = 0; i < src.imagew; i++)
for (j = 0; j < src.imageh; j++)
{
data += HCT.code[src.pDataAt(i)[j] + 1];
}
3、将编码后的图片数据写入文件
好了,编码已经完成了,下面需要把压缩后的图片数据写入文件了。首先我们把信息头、文件头、调色板的数据先写入:
//下面开始将信息头、文件头,这部分没有编码
FILE *f;
f = fopen("encode_Image.bmp","w+b");
if (f == NULL)
return false;
BITMAPFILEHEADER fh;
BITMAPINFOHEADER ih;
memset(&ih, 0, sizeof(BITMAPINFOHEADER));
fh.bfType = 0x4d42;
fh.bfReserved1 = fh.bfReserved2 = 0;
fh.bfOffBits = sizeof(BITMAPFILEHEADER) + sizeof(BITMAPINFOHEADER) + ((src.iYRGBnum == 1) ? 256 * sizeof(RGBQUAD) : 0);
ih.biSize = 40;
ih.biPlanes = 1;
ih.biWidth = src.imagew;
ih.biHeight = src.imageh;
ih.biBitCount = 8 * src.iYRGBnum;
int w4b = (src.imagew * src.iYRGBnum + 3) / 4 * 4;
ih.biSizeImage = ih.biHeight * w4b;
fh.bfSize = fh.bfOffBits + ih.biSizeImage;
fwrite(&fh, sizeof(BITMAPFILEHEADER), 1, f);
fwrite(&ih, sizeof(BITMAPINFOHEADER), 1, f);
if (src.iYRGBnum == 1) //因为我们只处理灰度图像,所以这里会需要写入调色板
fwrite(src.palette, sizeof(RGBQUAD), 256, f);
fclose(f);
下面,我们要处理像素数据域的数据了。刚才第二步中我们把用赫夫曼编码表示的像素数据存到了string data;
这个大字符串中,下面我们要把这个’1’、'0’构成的字符串,当成10二进制码写入文件了。
要完成这个操作还是有些困难的,我在尝试将字符转换为bitset
无果后,还是去搜寻了一番大神的博客,真的被我找到了一个可以实现我们目标的函数,这里放上博客地址,感谢这位大大。
需要注意的是,我们是在写完信息头、文件头等信息的文件中续写像素数据,因此需要将"wb"
模式改为"ab"
模式。
writeBit(data, "encode_Image.bmp"); //记得修改writeBit参数为string
可以用sublime text等软件查看其数据:(灰度图像有54字节的信息头和文件头,1000字节的调色板,因此像素数据域从第45行倒数第二个字节开始,如果你和我一样用的1.bmp这张图,其值应为4c)
至此,我们已经完成了对原始图像的赫夫曼编码。
4、读取压缩后的图像文件和赫夫曼编码表
之后是解码的操作,首先我们要从文件读回赫夫曼编码表,以及图像数据,再根据编码表还原出原始数据。
4.1 读取赫夫曼编码表
回顾一下我们写的编码表文件,其格式为一行之中,码字'\t'码长'\n'
,我们需要读取这个文件来恢复出一个HuffmanCodeTable HCT;
根据上述的需求,我们可以按行读取,将遇到'\t'
之前的部分赋值给码字,将'\t'
之后的部分,转换为数字,赋值给码长。代码如下:
bool readHuffmanTable(const char* cFilename,HuffmanCodeTable &HCT) //从文件读取赫夫曼码表
{
HCT.code = (HuffmanCode)malloc((256 + 1) * sizeof(char*)); //为HCT分配空间
HCT.length = (int*)malloc((256 + 1) * sizeof(int));
memset(HCT.length, 0, 257 * sizeof(int));
FILE *f;
f = fopen(cFilename, "r");
if (f == NULL)
return false;
char buffer[50];
int i=1, j=0, length;
while (fgets(buffer, 50, f) != NULL)
{
length = strlen(buffer);
HCT.code[i] = (char*)malloc(sizeof(char)*(length+1));
for (j = 0; j < length; j++)
{
if (buffer[j] == '\t')
{
HCT.code[i][j] = '\0';
while (buffer[++j] != '\n')
{
HCT.length[i] += (int)((buffer[j] - '0')*(pow((double)10, (double)(length - (j + 1) - 1))));
}
break;
}
else
HCT.code[i][j] = buffer[j];
}
i++;
}
return true;
}
可以注意到,对于码长的赋值,我自己推导了一个公式。举一个简单的例子,很容易搞清楚这个公式。
假设,某一行内容如下:1111100110 10,那么该行的长度length为14('\t'
和结尾的'\n'
不要忘记)。读到'\t'
时的下标j为10,说明下一位开始表示的就是码长了。但是码长不一定只有一位,因此我们写一个内部循环来求码长。
while (buffer[++j] != '\n')
j=11时指向的是1,这个1应该是个位、十位还是百位呢(如果有的话,实际上百位不太可能)?显然,我们可以通过length与j的值来推导出指数。j=11时指向的'1'
,表示的是:
同理,j=12时指向的'0'
表示的是:
因此我们就得到了这个通式:
+=
4.2 读取图像数据
就像写入时一样,我们将信息头、文件头、调色板先读取出来,确定图像的长、宽、位数等基本信息。
FILE *f;
if ((f = fopen(cFilename, "r+b")) == NULL)
{
return false;
}
//文件头、信息头的读取,和LoadBMPFILE函数一致
BITMAPFILEHEADER fh; //文件头
BITMAPINFOHEADER ih; //信息头
fread(&fh, sizeof(BITMAPFILEHEADER), 1, f); //从文件中读取文件头
if (fh.bfType != 0x4d42) //如果不是"BM"文件,直接退出
{
fclose(f);
return false;
}
fread(&ih, sizeof(BITMAPINFOHEADER), 1, f); //从文件中读取信息头
if ((ih.biBitCount != 8) && (ih.biBitCount != 24)) //如果不是8位灰度图,或24位彩图,则退出
{
fclose(f);
return false;
}
des.iYRGBnum = ih.biBitCount / 8; //1为灰度图,3为彩图
des.imagew = ih.biWidth;
des.imageh = ih.biHeight;
if (!des.AllocateMem()) //如果分配内存失败,退出
{
fclose(f);
return false;
}
if (des.iYRGBnum == 1) //灰度图需要调色板
fread(des.palette, sizeof(RGBQUAD), 256, f);
//信息头、文件头读取完毕
下面我们要处理像素数据了,之前写入时我们把字符串以二进制输出,现在我们要把这些二进制的01重新存为字符串。好在这学期的另一门课《计算机组成原理》中我曾经做过将float型按位输出以检验是否符合IEEE754标准的实验,因此我有一个函数可以将输入的东西以二进制输出。
char* Byte2Binary(unsigned char buffer) //将一个字节转换成8bit保存到字符串
{
char *s = (char*)malloc(sizeof(char) * 9); //返回字符串,记得给'\0'预留位置
s[0] = '\0';
unsigned char *p, ch;
int i, j;
p = (unsigned char *)(&buffer);
for (i = sizeof(buffer) - 1; i >= 0; i--)
{
ch = *(p + i);
for (j = 0; j < 8; j++)
{
if ((ch << j) & 0x80) //0x80即128,二进制表示1000 0000,这里操作目的是取最高位
s[j] = '1';
else
s[j] = '0';
}
s[j] = '\0';
}
return s;
}
用这个函数,我们只需要按一个字节的大小(恰好是一个unsigned char)来读取文件,再将其返回的字符数组拼接起来即可。
unsigned char buf; //注意char是1字节,char*是4字节
int rc;
string data; //获得编码数据
while ((rc = fread(&buf, sizeof(unsigned char), 1, f) != 0))
{
data += Byte2Binary(buf); //按字节读入,按位输出,存入字符串data
}
5、利用赫夫曼编码表还原图像文件
这一部分应该还不是最终版,因为我很不满意现在的执行效率,将样例图片还原需要一分钟的时间,这显然是不合格的。
简单说一下目前的思路,我们有了像素数据合并成的一个大字符串data,赫夫曼编码是不定长编码,每次取几个值来核对编码呢?这就需要我们求一下码表中的最短码长和最长码长。
int min = HCT.length[1], max = HCT.length[1]; //求最短码长和最大码长
for (i = 1; i <= 256; i++)
{
if (HCT.length[i] < min)
{
min = HCT.length[i];
}
if (HCT.length[i] > max)
{
max = HCT.length[i];
}
}
好了,现在我们从最短码长的码字开始匹配,一直匹配到最长码长的码字,如果有匹配的,那么这个像素值就解码成功了,如果全部走下来没有匹配的,说明编码数据有问题。
int index = 0; //字符串的下标标识
int flag;
for(i=0;i<des.imagew;i++) // 解码,给每个像素赋值
for (j = 0; j < des.imageh; j++)
{
flag = 0; //标志置0
if ((data.length() - (index+1)) >= min) //当剩余字符数大于最小码长时
{
for (k = min; k <= max; k++) //从最小码长到最大码长寻找对应码字
{
for (l = 1; l <= 256; l++)
{
if (HCT.length[l] == k) //对码长为k的码字,做比较
{
string temp;
temp += HCT.code[l];
if (!data.substr(index, k).compare(temp)) //注意string的compare方法,相等返回0
{
des.pDataAt(i)[j] = l;
flag = 1; //成功找到解码,标志置1
break;
}
}
}
if (flag == 1) //如果已经找到了像素值,就不再遍历
{
index += k; //下标移动
break;
}
}
if (flag == 0)
{
printf("编码数据有误!\n");
return false;
}
}
}
des.SaveBMPFILE("Recover.bmp");
这个方法简单直观,当然也很慢。显然码长并不是连续的,因此从最短码长遍历到最长码长做了很多无用功,我在想出优化方式后会对博客进行修改。
结语
这一次的实验终于是充满挑战性了,其实这篇博客已经把难点都规避掉了,因为真正的难点在于无从下手,在于走上了错误的解决道路却不自知,我和另外两名同学一起,花了整整三天时间,走遍了岔路才最终得到了结果。当然,收获确实是满满,就是头发又稀了不少。。。
本次实验完整的项目文件与代码可见我的gitee。