欢迎您访问程序员文章站本站旨在为大家提供分享程序员计算机编程知识!
您现在的位置是: 首页

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

程序员文章站 2022-07-14 21:55:36
...

一、实验原理

1、JPEG图像压缩标准基本介绍

      JPEG 是Joint Photographic Experts Group(联合图像专家小组)的缩写,是第一个国际图像压缩标准。JPEG图像压缩算法能够在提供良好的压缩性能的同时,具有比较好的重建质量,被广泛应用于图像、视频处理领域。

      JPEG本身只有描述如何将一个影像转换为字节的数据串流(streaming),但并没有说明这些字节如何在任何特定的储存媒体上被封存起来。.jpeg/.jpg是最常用的图像文件格式是一种有损压缩格式,能够将图像压缩在很小的储存空间,在获得极高的压缩率的同时能展现十分丰富生动的图像。

    JPEG格式压缩的主要是高频信息,对色彩的信息保留较好,适合应用于互联网,可减少图像的传输时间,可以支持24bit真彩色,也普遍应用于需要连续色调的图像。

2、JPEG压缩编码算法介绍

JPEG编码流程图如下:

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

编码步骤为:

(1)颜色空间转换、采样、零偏置

由于JPEG采用的是YUV颜色空间,首先将输入图片转换为YUV颜色空间,并按一定的采样格式进行采样。然后对Y分量进行电平偏移转换,可以对Y减128操作,使Y、U、V电平取值范围都在-128~127之间。

(2)分块

分别对Y、U、V三个分量的数据进行8*8分块,编码时,按从左至右,从上至下的顺序依次读取一个8*8块,对其进行DCT变换、量化、编码后,再读取下一个8*8块。

    若图像宽高不是8的整数倍,还需要通过复制相邻像素值补成8的整数倍。

(3)DCT变换

JPEG中采用DCT(离散余弦变换),DCT变换具有能量守恒、能量集中和去相关的特性。这样可以将能量集中在低频区域,对低频和高频各系数采用不同的量化比特数进行单独量化,不会损伤过多的效率。

(4)zig-zag扫描排序、量化

对每个8*8块的数据进行zig-zag扫描排序,重新排序后的8*8数组里的系数值按低频到高频依次排列。再对每个8*8系数块进行量化,由于人眼对亮度信号比对色差信号更敏感,所以采用亮度量化表和色差量化表对其分别处理。又由于人眼对低频比对高频敏感,所以对低频系数细量化,对高频粗量化。

zig-zag扫描:                                       两张量化表:

数据压缩实验五 JPEG原理分析及JPEG解码器的调试数据压缩实验五 JPEG原理分析及JPEG解码器的调试

(5)DC系数的差分脉冲调制编码

由于8*8 的图像块经过 DCT 变换之后得到的 DC 系数有两个特点:

(1)系数的数值比较大;

(2)相邻的 8*8 图像块的 DC 系数值变化不大:冗余;

所以对DC系数采用DPCM编码,即对相邻图像块之间量化DC系数的差值DIFF进行无损编码。

(6)AC系数的游程编码

因为DCT变换后AC系数按zig-zag排序后可以出现很多连零的机会,采用游程编码(RLC)可以进一步降低数据的传输量。

(7)熵编码

在得到 DC 系数的中间格式和 AC 系数的中间格式之后,为进一步压缩图象数据,有必要对两者进行熵编码。JPEG算法中通常采用Huffman编码方法。JPEG*采用了四张Huffman码表:亮度DC、亮度AC、色度DC、色度AC。

 

 3、JPEG文件格式介绍

       JPEG文件常用的格式称为JPEG文件交换格式(JPEG File Interchange Format,JFIF)。

       JPEG在文件中以 Segment的形式组织 ,它具有以下特点:

(1) 均以0xFF开始,后跟1 byte的 Marker 和 2 byte 的 Segment length(包含表示Length本身所占用的2 byte ,不含“0xFF” + “Marker” 占用的 2 byte );
(2) 采用 Motorola序(相对于Intel 序),即保存时高位在前,低位在后 ;

(3)Data 部分中, 0xFF后若为0x00,则跳过此字节不予处理。

 

       JPEG中常用的图像标记有:

(1)SOI ,Start of Image,   图像开始
  标记代码    2字节    固定值0xFFD8
(2)EOI, End of Image,   图像结束 

  标记代码    2字节    固定值0xFFD9

(3)APP0 应用程序保留标记0

  标记代码    2字节    固定值0xFFE0
包含9个具体字段:

 ① 数据长度                      2字节       ①~⑨9个字段的总长度

 ② 标识符                          5字节        固定值0x4A46494600,即字符串“JFIF0”

 ③ 版本号                          2字节        一般是0x0102,表示JFIF的版本号1.2  

 ④ X和Y的密度单位          1字节        只有三个值可选   0:无单位;1:点数/英寸;2:点数/厘米  

 ⑤ X方向像素密度             2字节        取值范围未知  

 ⑥ Y方向像素密度             2字节        取值范围未知  

 ⑦ 缩略图水平像素数目     1字节        取值范围未知  

 ⑧ 缩略图垂直像素数目     1字节        取值范围未知  

⑨ 缩略图RGB位图            长度可能是3的倍数     缩略图RGB位图数据
(4)DQT 定义量化表

标记代码     2字节    固定值0xFFDB
包含9个具体字段:

  ① 数据长度   2字节    字段①和多个字段②的总长度

  ② 量化表       数据长度-2字节

    a) 精度及量化表ID   1字节 
       高4位:精度,只有两个可选值    0:8位;1:16位    

       低4位:量化表ID,取值范围为0~3
    b) 表项         (64×(精度+1))字节
         例如8位精度的量化表,其表项长度为64×(0+1)=64字节
本标记段中,字段②可以重复出现,表示多个量化表,但最多只能出现4次

(5)SOF0   帧图像开始

标记代码      2字节     固定值0xFFC0
包含9个具体字段:  

 ① 数据长度    2字节       ①~⑥六个字段的总长度  

 ② 精度           1字节       每个数据样本的位数    通常是8位,一般软件都不支持 12位和16位  

 ③ 图像高度    2字节       图像高度(单位:像素)
 ④ 图像宽度    2字节      图像宽度(单位:像素)
 ⑤ 颜色分量数   1字节   只有3个数值可选   1:灰度图;3:YCrCb或YIQ;4:CMYK  

 而JFIF中使用YCrCb,故这里颜色分量数恒为3

 ⑥颜色分量信息   颜色分量数×3字节(通常为9字节)                                    

   a)颜色分量ID             1字节 
   b)水平/垂直采样因子 1字节 
      高4位:水平采样因子     低4位:垂直采样因子           

   c) 量化表  1字节      当前分量使用的量化表的ID

(6)DHT,定义哈夫曼表

标记代码      2字节      固定值0xFFC4
包含2个具体字段:
 ① 数据长度         2字节  
 ② huffman表       数据长度-2字节
    表ID和表类型   1字节  
    高4位:类型,只有两个值可选       0:DC直流;1:AC交流      

    低4位:哈夫曼表ID,

 注意,DC表和AC表分开编码

 不同位数的码字数量    16字节
编码内容 16个不同位数的码字数量之和(字节) 本标记段中,字段②可以重复出现(一般4次),也可以只出现1次。

(7)SOS      扫描开始

标记代码         2字节      固定值0xFFDA
包含2个具体字段:

 ①数据长度         2字节       ①~④两个字段的总长度

 ②颜色分量数      1字节      应该和SOF中的字段⑤的值相同,即:   1:灰度图是;3: YCrCb或YIQ;4:CMYK。
 ③颜色分量信息     

a) 颜色分量ID   1字节     

b) 直流/交流系数表号 1字节 
    高4位:直流分量使用的哈夫曼树编号         低4位:交流分量使用的哈夫曼树编号
 ④ 压缩图像数据        

 a)谱选择开始      1字节      固定值0x00        

 b)谱选择结束      1字节      固定值0x3F        

 c)谱选择             1字节      在基本JPEG中总为00


4、JPEG文件解码流程

JPEG文件解码过程与编码过程可逆,解码流程如下:

(1)读入文件的相关信息
(2)初步了解图像数据流的结构
(3)颜色分量单元的内部解码
(4)直流系数的差分编码
(5)反量化 & 反Zig-zag编码
(6)反离散余弦变换

(7)将YUV转换为需要的色彩空间并保存

二、实验流程及关键代码分析

1、JPEG算法系统中,所操作的对象数据具有明显的嵌套层次特点,我们可以采用分层的结构设计。

为了方便地对每个层次的数据进行操作,我们可以为每个层次数据设计一个结构体。

(1)JPEG解码时需要重建Huffman码表,我们可以定义一个Huffman码表结构体:

struct huffman_table
{
  short int lookup[HUFFMAN_HASH_SIZE];///建立查找表和权值的关系
  unsigned char code_size[HUFFMAN_HASH_SIZE];///建立权值和码长的关系
  uint16_t slowtable[16-HUFFMAN_HASH_NBITS][256];
};


(2)为了方便对每个8*8宏块进行Huffman解码,定义8*8块结构体:

struct component 
{
  unsigned int Hfactor;///水平采样因子
  unsigned int Vfactor;///垂直采样因子
  float *Q_table;	///指向量化表的指针
  struct huffman_table *AC_table;///AC码表
  struct huffman_table *DC_table;///DC码表
  short int previous_DC;	///前一个块的DC系数
  short int DCT[64];		///每个块的DCT系数
#if SANITY_CHECK
  unsigned int cid;
#endif
};


(3)再定义一个解码数据流结构体:

struct jdec_private
{
  uint8_t *components[COMPONENTS];///#define COMPONENTS	3 定义指针数组指向YUV三种块结构体
  unsigned int width, height;	///图像的宽、高
  unsigned int flags;

  const unsigned char *stream_begin, *stream_end;///解码位置标记指针
  unsigned int stream_length;

  const unsigned char *stream;	///当前解码流指针
  unsigned int reservoir, nbits_in_reservoir;

  struct component component_infos[COMPONENTS];///定义块结构体数组
  float Q_tables[COMPONENTS][64];		///定义YUV三张量化表
  struct huffman_table HTDC[HUFFMAN_TABLES];	///DC Huffman码表
  struct huffman_table HTAC[HUFFMAN_TABLES];	///AC Huffman码表  
  int default_huffman_table_initialized;
  int restart_interval;
  int restarts_to_go;				
  int last_rst_marker_seen;			

  uint8_t Y[64*4], Cr[64], Cb[64];///存放YUV分量值的数组

  jmp_buf jump_state;
  uint8_t *plane[COMPONENTS];
};

2、本实验中已经给出了JPEG解码的基本程序,我们可以采用单步调试的方式跟踪程序中数据的流向,掌握程序设计的整体框架。

读入文件名及输出格式后,从convert_one_image开始解码:

 

int convert_one_image(const char *infilename, const char *outfilename, int output_format)
{
  FILE *fp;
  unsigned int length_of_file;
  unsigned int width, height;
  unsigned char *buf;
  struct jdec_private *jdec;
  unsigned char *components[3];
 

  /* 将jpeg文件读入buf缓冲区 */
  fp = fopen(infilename, "rb");
  if (fp == NULL)
    exitmessage("Cannot open filename\n");
  length_of_file = filesize(fp);
  buf = (unsigned char *)malloc(length_of_file + 4);
  if (buf == NULL)
    exitmessage("Not enough memory for loading file\n");
  fread(buf, length_of_file, 1, fp);
  fclose(fp);

  jdec = tinyjpeg_init();///初始化
  if (jdec == NULL)
    exitmessage("Not enough memory to alloc the structure need for decompressing\n");

  if (tinyjpeg_parse_header(jdec, buf, length_of_file)<0)///解析头部信息
    exitmessage(tinyjpeg_get_errorstring(jdec));

  tinyjpeg_get_size(jdec, &width, &height);  ///获得图像宽高
 

  snprintf(error_string, sizeof(error_string),"Decoding JPEG image...\n");
  if (tinyjpeg_decode(jdec, output_format) < 0)///解码实际数据
    exitmessage(tinyjpeg_get_errorstring(jdec));

  tinyjpeg_get_components(jdec, components);

  /*保存为输出格式*/
  switch (output_format)
   {
    case TINYJPEG_FMT_RGB24:
    case TINYJPEG_FMT_BGR24:
      write_tga(outfilename, output_format, width, height, components);
      break;
    case TINYJPEG_FMT_YUV420P:
      write_yuv(outfilename, width, height, components);///输出成YUV420P文件
      break;
    case TINYJPEG_FMT_GREY:
      write_pgm(outfilename, width, height, components);///输出成灰度文件
	  break;
   }

  tinyjpeg_free(jdec);
  free(buf);
   
  return 0;
}


 

/*初始化函数*/
struct jdec_private *tinyjpeg_init(void)
{
  struct jdec_private *priv;
 
  priv = (struct jdec_private *)calloc(1, sizeof(struct jdec_private));
  if (priv == NULL)
    return NULL;
  return priv;
}


 

/*解析JPEG文件头*/
int tinyjpeg_parse_header(struct jdec_private *priv, const unsigned char *buf, unsigned int size)
{
  int ret;

  /* Identify the file */
  if ((buf[0] != 0xFF) || (buf[1] != SOI))//判断是否为JPEG文件
    snprintf(error_string, sizeof(error_string),"Not a JPG file ?\n");

  priv->stream_begin = buf+2;//begin指针前移2字节
  priv->stream_length = size-2;
  priv->stream_end = priv->stream_begin + priv->stream_length;
///开始解析JFIF
  ret = parse_JFIF(priv, priv->stream_begin);

  return ret;
}


 

/*解析各种不同的标签*/
static int parse_JFIF(struct jdec_private *priv, const unsigned char *stream)
{
  int chuck_len;
  int marker;
  int sos_marker_found = 0;
  int dht_marker_found = 0;
  const unsigned char *next_chunck;

  /* Parse marker */
  while (!sos_marker_found)
   {
     if (*stream++ != 0xff)
       goto bogus_jpeg_format;
     while (*stream == 0xff)
       stream++;

     marker = *stream++;///marker是跳过0xff字节后一个字节的内容
     chuck_len = be16_to_cpu(stream);///获取segment length
     next_chunck = stream + chuck_len;///指向下一个chunk的指针
     switch (marker)///判断是哪种marker标签
      {
       case SOF:
	 if (parse_SOF(priv, stream) < 0)
	   return -1;
	 break;
       case DQT:
	 if (parse_DQT(priv, stream) < 0)
	   return -1;
	 break;
       case SOS:
	 if (parse_SOS(priv, stream) < 0)
	   return -1;
	 sos_marker_found = 1;
	 break;
       case DHT:
	 if (parse_DHT(priv, stream) < 0)
	   return -1;
	 dht_marker_found = 1;
	 break;
       case DRI:
	 if (parse_DRI(priv, stream) < 0)
	   return -1;
	 break;
       default:
	 break;
      }
     stream = next_chunck;///解析下一个segment
   }
  。。。。。。}

解析DQT:

  • 得到量化表长度(可能包含多张)
  • 得到量化表的精度
  • 得到及检查量化表的序号(只能是0 —— 3)
  • 得到量化表内容( 64 个数据)
/*解析DQT*/
static int parse_DQT(struct jdec_private *priv, const unsigned char *stream)
{
  int qi;
  float *table;
  const unsigned char *dqt_block_end;
  dqt_block_end = stream + be16_to_cpu(stream);//得到量化表长度(可能包含多张量化表)
  stream += 2;	/* Skip length */

  while (stream < dqt_block_end)//检查是否还有表
   {
     qi = *stream++;
#if SANITY_CHECK
	 /*得到量化表的精度(高4位)*/
     if (qi>>4)
       snprintf(error_string, sizeof(error_string),"16 bits quantization table is not supported\n");
     /*得到量化表序号(低4位)*/
	 if (qi>4)
       snprintf(error_string, sizeof(error_string),"No more 4 quantization table is supported (got %d)\n", qi);
#endif
     table = priv->Q_tables[qi];
     build_quantization_table(table, stream);//得到量化表内容(64个数据)
     stream += 64;
   }
  return 0;
}


 


量化表以zig-zag排序存储:

/*zig-zag排序*/
static const unsigned char zigzag[64] = 
{
   0,  1,  5,  6, 14, 15, 27, 28,
   2,  4,  7, 13, 16, 26, 29, 42,
   3,  8, 12, 17, 25, 30, 41, 43,
   9, 11, 18, 24, 31, 40, 44, 53,
  10, 19, 23, 32, 39, 45, 52, 54,
  20, 22, 33, 38, 46, 51, 55, 60,
  21, 34, 37, 47, 50, 56, 59, 61,
  35, 36, 48, 49, 57, 58, 62, 63
};


 解析SOF:

  • 得到每个sample的比特数 、长宽、颜色分量数
  •  得到每个颜色分量的 ID、水平采样因子、垂直采样因子、使用的量化表序号(与 DQT中序号对应)
/*解析SOF*/
static int parse_SOF(struct jdec_private *priv, const unsigned char *stream)
{
  int i, width, height, nr_components, cid, sampling_factor;
  int Q_table;
  struct component *c;
  print_SOF(stream);

  height = be16_to_cpu(stream+3);//获取图像高度
  width  = be16_to_cpu(stream+5);//获取图像宽度
  nr_components = stream[7];//得到颜色分量数
  stream += 8;
  for (i=0; i<nr_components; i++) {//得到每个分量的信息
     cid = *stream++;//该分量的ID
     sampling_factor = *stream++;//该分量的水平、垂直采样因子
     Q_table = *stream++;//该分量使用的量化表序号
     c = &priv->component_infos[i];
     c->Vfactor = sampling_factor&0xf;//垂直采样因子(低4位)
     c->Hfactor = sampling_factor>>4;//水平采样因子(高4位)
     c->Q_table = priv->Q_tables[Q_table];

  }
  priv->width = width;
  priv->height = height;

  return 0;
}


解析DHT:

  • 得到 Huffman表的类型( AC 、DC )、序号
  • 依据数据重建Huffman表
/*解析DHT*/
static int parse_DHT(struct jdec_private *priv, const unsigned char *stream)
{
  unsigned int count, i;
  unsigned char huff_bits[17];
  int length, index;

  length = be16_to_cpu(stream) - 2;///表长(可能含有多张表)
  stream += 2;	/* Skip length */

  while (length>0) {//是否还有表
     index = *stream++;

     /* We need to calculate the number of bytes 'vals' will takes */
     huff_bits[0] = 0;
     count = 0;
     for (i=1; i<17; i++) {
	huff_bits[i] = *stream++;
	count += huff_bits[i];
     }
     if (index & 0xf0 )//AC表
       build_huffman_table(huff_bits, stream, &priv->HTAC[index&0xf]);
     else//DC表
       build_huffman_table(huff_bits, stream, &priv->HTDC[index&0xf]);

     length -= 1;
     length -= 16;
     length -= count;
     stream += count;
  }
  return 0;
}


解析SOS:

得到解析每个颜色分量的 DC 、AC 值所使用的Huffman 表序号(与 DHT中序号对应)

/*解析SOS*/
static int parse_SOS(struct jdec_private *priv, const unsigned char *stream)
{
  unsigned int i, cid, table;
  unsigned int nr_components = stream[2];//得到颜色分量数
  stream += 3;
  for (i=0;i<nr_components;i++) {//对所有分量
     cid = *stream++;//得到ID
     table = *stream++;//Huffman表序号,高四位:DC;低四位:AC
     priv->component_infos[i].AC_table = &priv->HTAC[table&0xf];
     priv->component_infos[i].DC_table = &priv->HTDC[table>>4];
  }
  priv->stream = stream+3;
  return 0;
}



依据每个分量的水平、垂直采样因子MCU的大小,并得到每个MCU中 8*8块的个数:

 xstride_by_mcu = ystride_by_mcu = 8;///初始化每个MCU的宽、高为8,即采样格式4;4;4
  if ((priv->component_infos[cY].Hfactor | priv->component_infos[cY].Vfactor) == 1) {//水平、垂直采样因子都为1
     decode_MCU = decode_mcu_table[0];//每个MCU中1*1个Y分量块
     convert_to_pixfmt = colorspace_array_conv[0];
#if TRACE
     fprintf(p_trace,"Use decode 1x1 sampling\n");
	 fflush(p_trace);
#endif
  } else if (priv->component_infos[cY].Hfactor == 1) {//水平采样因子1,垂直采样因子2
     decode_MCU = decode_mcu_table[1];//每个MCU1*2个Y分量块
     convert_to_pixfmt = colorspace_array_conv[1];
     ystride_by_mcu = 16;//每个MCU的高为16
#if TRACE
     fprintf(p_trace,"Use decode 1x2 sampling (not supported)\n");
	 fflush(p_trace);
#endif
  } else if (priv->component_infos[cY].Vfactor == 2) {//水平采样因子2,垂直采样因子2
     decode_MCU = decode_mcu_table[3];//每个MCU2*2个Y分量块
     convert_to_pixfmt = colorspace_array_conv[3];
     xstride_by_mcu = 16;//每个MCU的宽为16
     ystride_by_mcu = 16;//每个MCU的高为16
#if TRACE 
	 fprintf(p_trace,"Use decode 2x2 sampling\n");
	 fflush(p_trace);
#endif
  } else {
     decode_MCU = decode_mcu_table[2];//水平采样因子2,垂直采样因子1,每个MCU2*1个Y分量块
     convert_to_pixfmt = colorspace_array_conv[2];
     xstride_by_mcu = 16;//每个MCU的宽为16
#if TRACE
     fprintf(p_trace,"Use decode 2x1 sampling\n");
	 fflush(p_trace);
#endif
  }



对每个MCU解码(依照各分量水平、垂直采样因子对 MCU中每个分量宏块解码)
1 对每个宏块进行 Huffman解码,得到DCT系数
2 对每个宏块的 DCT系数进行 IDCT,得到  Y、Cb 、Cr
3 遇到 Segment Marker RST时,清空之前的DC  DCT系数。

/*
 * Decode all the 3 components for 1x1 
 */
static void decode_MCU_1x1_3planes(struct jdec_private *priv)//4:4:4
{
  // Y
  process_Huffman_data_unit(priv, cY);//以8*8块为单位进行Huffman解码
  IDCT(&priv->component_infos[cY], priv->Y, 8);//对得到的DCT系数进行IDCT
  
  // Cb
  process_Huffman_data_unit(priv, cCb);
  IDCT(&priv->component_infos[cCb], priv->Cb, 8);

  // Cr
  process_Huffman_data_unit(priv, cCr);
  IDCT(&priv->component_infos[cCr], priv->Cr, 8);
}
/*
 * Decode a 2x1
 *  .-------.
 *  | 1 | 2 |
 *  `-------'
 */
static void decode_MCU_2x1_3planes(struct jdec_private *priv)//4:2:2
{
  // Y
  process_Huffman_data_unit(priv, cY);//第1个块
  IDCT(&priv->component_infos[cY], priv->Y, 16);
  process_Huffman_data_unit(priv, cY);//第2个块
  IDCT(&priv->component_infos[cY], priv->Y+8, 16);

  // Cb
  process_Huffman_data_unit(priv, cCb);
  IDCT(&priv->component_infos[cCb], priv->Cb, 8);

  // Cr
  process_Huffman_data_unit(priv, cCr);
  IDCT(&priv->component_infos[cCr], priv->Cr, 8);
}
/*
 * Decode a 2x2
 *  .-------.
 *  | 1 | 2 |
 *  |---+---|
 *  | 3 | 4 |
 *  `-------'
 */
static void decode_MCU_2x2_3planes(struct jdec_private *priv)//4:2:0
{
  // Y
	///////////////////////////////////
  process_Huffman_data_unit(priv, cY);//第1个块
  IDCT(&priv->component_infos[cY], priv->Y, 16);
  process_Huffman_data_unit(priv, cY);//第2个块
  IDCT(&priv->component_infos[cY], priv->Y+8, 16);
  process_Huffman_data_unit(priv, cY);//第3个块
  IDCT(&priv->component_infos[cY], priv->Y+64*2, 16);
  process_Huffman_data_unit(priv, cY);//第4个块
  IDCT(&priv->component_infos[cY], priv->Y+64*2+8, 16);

  // Cb
  process_Huffman_data_unit(priv, cCb);
  IDCT(&priv->component_infos[cCb], priv->Cb, 8);

  // Cr
  process_Huffman_data_unit(priv, cCr);
  IDCT(&priv->component_infos[cCr], priv->Cr, 8);
}
/*
 * Decode a 1x2 mcu
 *  .---.
 *  | 1 |
 *  |---|
 *  | 2 |
 *  `---'
 */
static void decode_MCU_1x2_3planes(struct jdec_private *priv)
{
  // Y
  process_Huffman_data_unit(priv, cY);
  IDCT(&priv->component_infos[cY], priv->Y, 8);
  process_Huffman_data_unit(priv, cY);
  IDCT(&priv->component_infos[cY], priv->Y+64, 8);

  // Cb
  process_Huffman_data_unit(priv, cCb);
  IDCT(&priv->component_infos[cCb], priv->Cb, 8);

  // Cr
  process_Huffman_data_unit(priv, cCr);
  IDCT(&priv->component_infos[cCr], priv->Cr, 8);
}


以8*8宏块为单位进行Huffman解码:

static void process_Huffman_data_unit(struct jdec_private *priv, int component)
 {unsigned char j;
  unsigned int huff_code;
  unsigned char size_val, count_0;

  struct component *c = &priv->component_infos[component];
  short int DCT[64];


  /* Initialize the DCT coef table */
  memset(DCT, 0, sizeof(DCT));

  /* DC coefficient decoding */
  huff_code = get_next_huffman_code(priv, c->DC_table);
  //trace("+ %x\n", huff_code);
  if (huff_code) {
     get_nbits(priv->reservoir, priv->nbits_in_reservoir, priv->stream, huff_code, DCT[0]);
	 //解码得到当前块与前一块DC系数的差值
     DCT[0] += c->previous_DC;//差值与前一DC系数相加,得到当前块DC系数
     c->previous_DC = DCT[0];//更新前一DC系数值
  } else {
     DCT[0] = c->previous_DC;
  }

  /* AC coefficient decoding */
  j = 1;
  while (j<64)
   {
     huff_code = get_next_huffman_code(priv, c->AC_table);
     //trace("- %x\n", huff_code);

     size_val = huff_code & 0xF;//幅度
     count_0 = huff_code >> 4;//零游程长度

     if (size_val == 0)//0不是一个有效的Amplitude值,这里做零游程标志
      { /* RLE */
	if (count_0 == 0)
	  break;	/* EOB found, go out */
	else if (count_0 == 0xF)
	  j += 16;	/* skip 16 zeros */
      }
     else
      {
	j += count_0;	/* skip count_0 zeroes */
	if (__unlikely(j >= 64))//出错了
	 {
	   snprintf(error_string, sizeof(error_string), "Bad huffman data (buffer overflow)");
	   break;
	 }
	get_nbits(priv->reservoir, priv->nbits_in_reservoir, priv->stream, size_val, DCT[j]);
	//查表得到AC DCT系数
	j++;
      }
   }

  for (j = 0; j < 64; j++)
    c->DCT[j] = DCT[zigzag[j]];//以zig-zag序保存

}


循环解完所有MCU,解码结束

 for (y=0; y < priv->height/ystride_by_mcu; y++)//行循环
   {
     //trace("Decoding row %d\n", y);
     priv->plane[0] = priv->components[0] + (y * bytes_per_blocklines[0]);
     priv->plane[1] = priv->components[1] + (y * bytes_per_blocklines[1]);
     priv->plane[2] = priv->components[2] + (y * bytes_per_blocklines[2]);
     for (x=0; x < priv->width; x+=xstride_by_mcu)//列循环
      {
	decode_MCU(priv);//对一个MCU解码(Huffman解码+IDCT)

	convert_to_pixfmt(priv);
	priv->plane[0] += bytes_per_mcu[0];
	priv->plane[1] += bytes_per_mcu[1];
	priv->plane[2] += bytes_per_mcu[2];
	if (priv->restarts_to_go>0)
	 {
	   priv->restarts_to_go--;
	   if (priv->restarts_to_go == 0)
	    {
	      priv->stream -= (priv->nbits_in_reservoir/8);
	      resync(priv);//清空preDC(所有颜色分量)
	      if (find_next_rst_marker(priv) < 0)//查找RST标记
		return -1;
	    }
	 }
      }
   }


三、实验调试与结果分析

1、将输出文件保存为可供YUVviewer观看的YUV文件:

static void write_yuv(const char *filename, int width, int height, unsigned char **components)
{
  FILE *F;
  char temp[1024];
  snprintf(temp, 1024, "%s.yuv", filename);//输出一个.yuv文件
  F = fopen(temp, "wb");//打开输出文件
  fwrite(components[0], width*height, 1,F);//将Y分量数据写入输出文件
  fwrite(components[1], width*height/4, 1,F);//将U分量数据写入输出文件
  fwrite(components[2], width*height/4, 1,F);//将V分量数据写入输出文件
  fclose(F);}

输入一张960*540的图片,输出yuv图像如下:

数据压缩实验五 JPEG原理分析及JPEG解码器的调试
由于图像高度为540,不是8的整数倍,解码时将高度压缩为536,这样就可以做8*8分块处理。所以输出的yuv图像下方出现了4行空值。

2、以txt文件输出所有的量化矩阵

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

定位到两个FF DB,说明有两张量化表;

红色框为量化表长度:00 43H表示67字节;

黄色框00:低位0表示量化表ID索引值,高位为量化精度,0表示8bit,即用一个字节表示一个量化系数。

///change by zgy for output DQT
void output_DQT(const char* infilename,FILE *DQT_txt)
{
FILE *fp;
int i,j;
int length_of_file;
unsigned char *buf;
fp = fopen(infilename, "rb");//打开输入文件
  length_of_file = filesize(fp);
  buf = (unsigned char *)malloc(length_of_file );
  fread(buf, length_of_file, 1, fp);//将输入文件读入缓冲区
  fclose(fp);

  for(i=0;i<length_of_file;i++)
  {
if(buf[i]==0xFF&&buf[i+1]==0xDB)//定位到DQT marker
{
	unsigned char*DQTbuf;
	DQTbuf=(unsigned char *)malloc(buf[i+2]*16*16+buf[i+3]-3);//量化表数据长度为67-3=64字节
	fprintf(DQT_txt,"DQT marker (length=%d)\n",buf[i+2]*16*16+buf[i+3]-3);
	for(j=0;j<buf[i+2]*16*16+buf[i+3]-3;j++)
	{
		DQTbuf[j]=buf[zigzag[j]+i+5];//量化表在JPEG文件中按zig-zag排序,查表变换后按从左至右正常序输出
		fprintf(DQT_txt,"%2.2x ",DQTbuf[j]);
		if((j+1)%8==0)
			fprintf(DQT_txt,"\n");
	}
	free(DQTbuf);
  }
  }
free(buf);
}
///end change

输出txt文件如下:

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

3、以txt文件输出所有的Huffman码表

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

定位到4个FF C4,表示有4张Huffman码表;

红色框的两字节表示码表长度:例如第一张码表00 1F为31字节;

黄色框的一字节表示码表ID和类型:分别按DC0、AC0、DC1、AC1顺序排列;

蓝色框的16个字节表示不同位数的码字数量,每个字节代表从长度1到长度16的码字个数。可以从DC0表中得出,一共有1+5+1*6=12个码字;

绿色框即为这12个码字对应的权值,共12字节。

 

在程序调试时,我们可以通过添加条件编译命令#if……#endif,输出一些中间结果,例如在这里我们可以用#if TRACE……#endif输出Huffman码表查看:

#define TRACE 1
#define  TRACEFILE "trace_jpeg.txt"

将TRACE变量定义成常量1,可以执行所有的条件编译指令,编译程序中所有满足条件的程序段;如果不再输出中间结果,只用将TRACE改成0即可。这样便于程序的调试、修改,便于实现模块化层次设计。

输出一些表头信息:

#if TRACE
     fprintf(p_trace,"Huffman table %s[%d] length=%d\n", (index&0xf0)?"AC":"DC", index&0xf, count);
	 fflush(p_trace);
#endif


在build_Huffman_table函数中添加条件编译命令输出重建的码表:

#if TRACE
     fprintf(p_trace,"val=%2.2x code=%8.8x codesize=%2.2d\n", val, code, code_size);
	 fflush(p_trace);
    #endif


输出的TRACEFILE文件如下:

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

3、输出DC图像并Huffman统计其概率分布

为了输出DC图像,应该找到每个8*8块解码后的DCT[0]即DC系数,依次写入DCbuf缓冲区,再输出DCimage文件。

我们可以在tinyjpeg_decode解码函数中添加代码输出DCimage文件:

///change by zgy for DCimage
 int i=0,j=0,n=0;
 int k=0;
 FILE *DC_fp;
  unsigned char *DCbuf;
  DC_fp=fopen("DCimage.yuv","wb");

   DCbuf=(unsigned char *)malloc(priv->height*priv->width/64);
  ///end change


由于每个块的大小为8*8=64个像素,而每个块DCT变换后只有一个DC系数值,即DCT[0],所以为DCbuf缓冲区开辟空间priv->height*priv->width/64。

 

按照解码流程我们可以很快找到每个块的DC系数在process_Huffman_data_unit函数中,通过DC差分系数的解码得到差分校正分量,与previous_DC相加得到当前DC系数值,再反量化,反zig-zag存储,赋值给DCT[0]。

问题1:怎样将这个DCT[0]值通过函数调用一层层传给DCbuf呢?

这里笔者起初犯了一个低级错误,直接在process_Huffman_data_unit中添加参数int DC,将DCT[0]赋给DC,又在decode_MCU函数中同样添加参数DC,再一层层将DC传递到tinyjpeg_decode函数中给DCbuf。这样显然是不正确的,函数中实参和形参的存储空间是相互独立的,形参值的改变不能影响实参。

这里就可以发挥结构体的作用了,在这一层层函数的调用中,函数的参数值都包含结构体指针priv,通过指针指向结构体,我们就可以改变priv结构体中的值。对此,我们可以在定义的jdec_private结构体中添加两个变量用来传递DC系数:

///change by zgy for DCimage
  int DCC;///定义变量传递每个8*8块的直流系数
  int DC[4];///定义数组传递2*2宏块格式的MCU的4个直流系数
  ///end change

process_Huffman_data_unit函数中将DCT[0]值传给DCC:

///change by zgy for DCimage
 priv->DCC=c->DCT[0];///将当前宏块的直流系数传给DCC
 ///end change


先考虑MCU Y分量为2*2格式的情况,调用一次decode_MCU函数可以得到4个块的DC值,所以再定义一个DC[4]的数组用来传递这4个DC系数:

static void decode_MCU_2x2_3planes(struct jdec_private *priv)
{
  // Y
  process_Huffman_data_unit(priv, cY); 
  priv->DC[0]=priv->DCC; ///change by zgy for DCimage
  IDCT(&priv->component_infos[cY], priv->Y, 16);
  process_Huffman_data_unit(priv, cY); 
   priv->DC[1]=priv->DCC; ///change by zgy for DCimage
  IDCT(&priv->component_infos[cY], priv->Y+8, 16);
  process_Huffman_data_unit(priv, cY); 
   priv->DC[2]=priv->DCC; ///change by zgy for DCimage
  IDCT(&priv->component_infos[cY], priv->Y+64*2, 16);
  process_Huffman_data_unit(priv, cY); 
   priv->DC[3]=priv->DCC; ///change by zgy for DCimage
  IDCT(&priv->component_infos[cY], priv->Y+64*2+8, 16);

  // Cb
  process_Huffman_data_unit(priv, cCb); 
  IDCT(&priv->component_infos[cCb], priv->Cb, 8);

  // Cr
  process_Huffman_data_unit(priv, cCr); 
  IDCT(&priv->component_infos[cCr], priv->Cr, 8);
}

得到一个MCU的4个DC值后,我们就可以将其写入DCbuf了。

问题2:一个MCU的DC值是按田字格顺序排列的,而DCbuf中数据应该按从左至右,从上往下排列,如何调整数据输出顺序呢?

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

我们可以另外再开辟两个缓冲区分别用来存放奇数行和偶数行的数据:

  unsigned char *DCbuf1;///偶数行直流系数缓冲区
  unsigned char *DCbuf2;///奇数行直流系数缓冲区
 DCbuf1=(unsigned char *)malloc(priv->width/8);
  DCbuf2=(unsigned char *)malloc(priv->width/8);

将每个MCU中的DC[0]、DC[1]写入奇数行缓冲区,DC[2]、DC[3]写入偶数行缓冲区,再将2个缓冲区的数据合并。

	DCbuf1[2*i]=priv->DC[0]/8;
	DCbuf1[2*i+1]=priv->DC[1]/8;
	DCbuf2[2*j]=priv->DC[2]/8;
	DCbuf2[2*j+1]=priv->DC[3]/8;
	i++;
	j++;
for(n=0;n<priv->width/8;n++)
	 {
		 DCbuf[k*priv->width/8+n]=DCbuf1[n];
		 DCbuf[(k+1)*priv->width/8+n]=DCbuf2[n];
	 }
	k+=2;
	i=0;
	j=0;


值得注意的是DC系数的取值范围,原图中Y分量的像素值范围是[0,255],DCT变换后,根据能量守恒特性,DC系数的取值范围应该在[0,8*255],为了准确输出DC图像,可以将DC值都除以8,使DC图像灰度值在[0,255]范围内。

输入一张2*2格式的图片,输出DCimage如下:

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

图像已有大致轮廓,但是像素值明显不对。将原图Y分量文件与DCimage输入Huffman编码器,对比概率分布图:

 

数据压缩实验五 JPEG原理分析及JPEG解码器的调试数据压缩实验五 JPEG原理分析及JPEG解码器的调试

DCimage灰度值集中在0附近,与原图概率分布不同,调试代码发现previous_DC初始为0,将其改成128*8,再输出图像。

static void resync(struct jdec_private *priv)
{
  int i;

  /* Init DC coefficients */
  for (i=0; i<COMPONENTS; i++)
     priv->component_infos[i].previous_DC = 128*8;///change by zgy for DCimage

  priv->reservoir = 0;
  priv->nbits_in_reservoir = 0;
  if (priv->restart_interval > 0)
    priv->restarts_to_go = priv->restart_interval;
  else
    priv->restarts_to_go = -1;
}

数据压缩实验五 JPEG原理分析及JPEG解码器的调试数据压缩实验五 JPEG原理分析及JPEG解码器的调试

可以看到输出的DCimage轮廓清晰,灰度分明,大致为原图的缩略图,其概率分布趋势与原图相同,只是因为previous_DC初始为128*8,除以8之后为128,使得灰度值压缩在128附近。


 问题3:对不同的MCU格式怎样处理?

static const decode_MCU_fct decode_mcu_3comp_table[4] = {///MCU宏块格式数组
   decode_MCU_1x1_3planes,
   decode_MCU_1x2_3planes,
   decode_MCU_2x1_3planes,
   decode_MCU_2x2_3planes,
};


针对不同采样格式的图片,一个MCU中所包含Y分量块的个数不同,每次decode_MCU后得到的DC系数个数也不同,需要分情况处理。

///change by zgy for DCimage
	if(decode_MCU==decode_mcu_table[0])
	{
		DCbuf[i]=priv->DC[0]/8;
		i++;
	}
	else if(decode_MCU==decode_mcu_table[1])
	{
		DCbuf1[i]=priv->DC[0]/8;
		DCbuf2[j]=priv->DC[2]/8;
		i++;
		j++;
		}
		else if(decode_MCU==decode_mcu_table[2])
		{
			DCbuf[i*2]=priv->DC[0]/8;
			DCbuf[i*2+1]=priv->DC[1]/8;
			i++;
		}
	else
		{
			DCbuf1[2*i]=priv->DC[0]/8;
			DCbuf1[2*i+1]=priv->DC[1]/8;
			DCbuf2[2*j]=priv->DC[2]/8;
			DCbuf2[2*j+1]=priv->DC[3]/8;
			i++;
			j++;
		}
 ///end change

其中1x2和2x2格式需要调整数据写入DCbuf的顺序:

///change by zgy for DCimage
	 if(decode_MCU==decode_mcu_table[1]||decode_MCU==decode_mcu_table[3])
	 {
	 for(n=0;n<priv->width/8;n++)
	 {
		 DCbuf[k*priv->width/8+n]=DCbuf1[n];
		 DCbuf[(k+1)*priv->width/8+n]=DCbuf2[n];
	 }
	k+=2;
	i=0;
	j=0;
	 }
		 ///end change


用不同采样格式图片测试:

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

数据压缩实验五 JPEG原理分析及JPEG解码器的调试

4、输出AC图像并Huffman统计其概率分布

对AC系数的操作与DC系数类似,在结构体中jdec_private中添加两个变量传递AC系数:

///change by zgy for ACimage
int ACC;
int AC[4];
  ///end change

在process_Huffman_data_unit函数中添加语句:

  ///change by zgy for ACimage
priv->ACC=c->DCT[1]+128;
  ///end change


由于AC系数表示高频分量,概率分布集中在0附近,加128之后便于输出。AC图像与概率分布图如下:

DCT[1]:

数据压缩实验五 JPEG原理分析及JPEG解码器的调试数据压缩实验五 JPEG原理分析及JPEG解码器的调试

DCT[9]:

数据压缩实验五 JPEG原理分析及JPEG解码器的调试数据压缩实验五 JPEG原理分析及JPEG解码器的调试

可以验证,DCT系数频率分量越高,概率分布越趋近于0.

5.高频截止对图像的影响

锐利的高频截止,会造成时域分散,也就是图像会产生“蚊子噪声”。而图像的压缩比越大,高频截止的成分越多。下面对不同图像质量的图像进行对比分析:

 

图像质量等级 图像 图像细节
数据压缩实验五 JPEG原理分析及JPEG解码器的调试 数据压缩实验五 JPEG原理分析及JPEG解码器的调试
中等 数据压缩实验五 JPEG原理分析及JPEG解码器的调试 数据压缩实验五 JPEG原理分析及JPEG解码器的调试
数据压缩实验五 JPEG原理分析及JPEG解码器的调试 数据压缩实验五 JPEG原理分析及JPEG解码器的调试

可以看到,图像质量越低,蚊子噪声越多,即时域分散现象越明显。

四、实验结论

通过对JPEG解码程序的调试、分析,更加深了我对一个算法型系统的理解。如何对一个完整的系统进行模块划分与结构化层次设计,可以帮助我们理清数据的流向,快速地获取到所需要的数据信息。通过适当添加条件编译命令输出中间结果,便于对程序的调试、修改。分析问题时要理清逻辑思路,通过分步调试跟踪数据的流向,找到需要优化的地方,将问题分解逐个击破。JPEG解码程序中还有许多细节地方值得我们借鉴学习。