图像编码是图像应用的基础，随着云时代和大数据时代的到来，网络中规模巨大的视频和图像数据为存储与数据传输带来了极大压力。多年来，对图像编码方法的研究从未停止，研究者始终致力于研究兼顾高图像质量与高压缩比的算法^[1-5]。传统的图像编码算法大多是利用图像自身像素之间的相关性进行编码，图像编码的国际标准有JPEG、JPEG 2000^[6]和HEVC/H.265^[7]。这些算法更多的适用于独立图像的编码，在云端拥有海量数据资源的情况下，极有可能找到与待编码图像具有相似部分的图像，传统的编码算法无法利用这种相关性进行编码，因此会产生很多不必要的冗余。如果把云端的海量数据资源看作一个大的图像库，利用库中已有图像与待编码图像具有相似特征的部分对图像进行编码，则能够大量节省编码所需的码流，达到图像压缩的目的。

利用相似的图像对原图像进行编码，首先要在库中找到这些相似图像。但即使图像库中的数据量再大，也很难找到整幅图像都与原图非常相似的图像，更多的仅是图像局部存在相同或相似的图像块，因此考虑用局部特征来对图像进行分块描述和匹配。David提出SIFT算法^[8]是一种基于尺度空间的图像局部特征提取算法，可以对局部图像中的极值点进行检测并描述。尺度不变性和旋转不变性使得SIFT算法对图像的缩放和旋转具有良好的适应性，特征点在尺度不同或方向不同的情况下也能准确匹配。近年来，SIFT算法也在图像识别^[9-11]、图像分类^[12-13]和图像拼接^[14-15]等方面取得了许多成果，这些都为基于SIFT算法的图像编码奠定了理论基础。

利用SIFT特征对图像进行编码，如果直接在空域进行特征匹配与重构，会使得重构质量很不理想^[16]。本文提出一种小波变换下SIFT特征匹配图像编码方法，通过小波变换^[17]，将图像分解为低频子带和高频子带。低频子带携带原图像大部分的颜色和结构信息，经过编码传到解码端，以保证重构图像的质量；高频子带携带细节和纹理信息，对高频部分提取特征，将其由像素值表达转为特征描述表达。在解码端根据特征描述，在库中找到匹配图像块完成高频子带的复原，再通过小波逆变换得到重构图像。在高频部分提取到的特征点数量远远小于直接在原图上提取特征点的数量，压缩低频子带占用的空间也远小于压缩原图所需的空间，因此本文算法可以在保证图像重构质量的同时达到很高的压缩比。

1 编码实现

图 1为编码实现的流程框图。通过小波变换将输入图像的低频部分和高频部分分解出来。编码系统也由此分为低频子带编码和高频子带编码两大部分。

1.1 小波变换

输入图像经过1次小波变换以后可以分解为1个低频子带(LL)和3个高频子带(LH, HL, HH)，其中，LL携带了原图中物体的结构信息和大部分的色彩信息，LH、HL和HH分别携带了水平、垂直和对角线方向的纹理信息。经过小波分解后的低频子带可以进一步被分解，第n层小波分解得到的子带可以表示为LL_n、LH_n、HL_n、HH_n。图 2是2层小波分解的示意图，可以看出，小波变换的层数越多，低频子带所占用的空间就越小，但同时低频子带携带的信息也越少。

本文采用的小波基为Haar小波，可以表示为

$ \psi \left( t \right) = \left\{ \begin{array}{l} 1,\;\;\;\;t \in \left[ {0,1/2} \right)\\ - 1,\;\;\;\;t \in \left[ {1/2,1} \right] \end{array} \right. $

(1)

由于所有实验用的输入图像均为3通道的彩色图像，而小波变换适用于单通道图像，因此需要分别对3个通道进行小波变换，得到每个通道的小波分解系数。即经过小波分解后，分别得到R、G、B三个通道的小波系数。

引入小波变换目的有2个：1)原图中的结构信息是全局信息，而纹理信息是局部信息，通过小波变换分解出高频部分的纹理信息，才能更好的应用SIFT这样的局部特征提取算法；2)小波分解将原图的信息重新分布，低频子带携带的信息较多，高频子带携带的信息较少。对低频子带采用低失真的编码算法，而高频子带采用高压缩比的编码算法，可以有效的提高编码效率。

1.2 低频子带编码

根据上文所述，低频子带要采用低失真的图像编码算法，经典的图像编码算法可用于低频子带的压缩。本文主要使用的是JPEG 2000对低频子带进行编码。

1.3 高频子带编码

经过n层小波变换后，共产生3×n个高频子带，为了得到全面的高频信息，首先要进行n次无低频部分的小波逆变换(将低频子带的值都取0，再结合高频子带的小波逆变换)，使得所有高频信息都汇集到一幅图像中，即高频图像。逆变换后得到的高频图像如图 3所示。

在获得高频图像后，对该图像提取SIFT特征点，通过特征提取过程可以找到图像中存在的边缘、角点及一些纹理区域，这些区域需要被保留。而对于没有特征的区域，通常是图像中的一些平坦区，如墙面，平静的水面等，这些区域则不会被保留。换句话说，特征提取这一过程选择性的保留了高频图像中的有用区域，从而有效的减少了冗余。经过这一过程，提取到的特征点可以表示为

$ K = \left\{ {{k_i}\left( {{x_i},{y_i},{s_i},{o_i}} \right),i = 1,2, \cdots ,n} \right\} $

(2)

式中：(x_i, y_i)是特征点的位置；s_i为特征点的尺度，表示是在尺度空间的哪一层提取到的特征点；o_i为特征点的主方向；n为提取到的特征点的数量。

特征点的位置信息在解码端将用于高频子带的复原，因此需要对其进行无损压缩，本文采用游程长度编码来编码位置信息，首先建立一个与高频图像尺度相同的矩阵P，并将矩阵中所有元素置0。然后根据特征点的位置坐标，将矩阵中的对应元素置1。即

$ \mathit{\boldsymbol{P}}\left( {{x_i},{y_i}} \right) = 1,i = 1,2, \cdots ,n $

(3)

经过这一步骤，可以得到一个二值矩阵，对这个二值矩阵进行二进制游程长度编码，得到的码流即为位置信息编码码流。

特征点的尺度和方向信息与位置信息一起用于计算特征描述子，SIFT算法的特征描述子是根据特征点周围区域像素之间的梯度关系计算出的一个128维的向量^[8]。不同的特征点之间可以通过特征描述子来进行匹配，最常用的方法是计算两个描述子之间的欧式距离：

$ {\rm{dis}}{{\rm{t}}_{ij}} = \sum\limits_{m = 1}^{128} {\sqrt {{{\left( {{d_{im}} - {d_{jm}}} \right)}^2}} } $

(4)

式中：d_i、d_j分别为两个特征点k_i和k_j的特征描述子，d_im表示描述子d_i中第m维的数值。

每张图上可以提取到成百上千个特征点，每个SIFT特征点的描述子都是一个128维的向量，如果直接对描述子本身的数据进行编码，会占用大量的比特数。因此，转而利用现有的图像库，将图像库中所有的特征点作为训练集，训练出一个规定大小的码书。码书的训练过程表示为

$ {C^D} = Q{\left( {{D^D}} \right)_{size\;c}} $

(5)

式中：C^D是由D^D经过矢量训练得到的码书，本文采用矢量量化^[18]的方法对所有高频子带提取到的特征描述子进行训练；D^D是图像库中所有特征点描述子的集合；size c是规定的码书大小。由此可见，C^D也是一个描述子的集合，这个集合中共有size c个特征描述子，在编码端，用高频图像提取到的特征点与码书中的特征点进行匹配，可以得到匹配的索引序号，将这个索引序号传到解码端，就可以在解码端利用同样的码书解码出特征描述子。

但是，由于描述子的特异性，如果单纯用码书数据来代替描述子，会使匹配产生较大的误差，因此除了对码书索引进行编码以外，还需根据式(6)计算描述子残差：

$ {\delta _{im}} = \left\{ \begin{array}{l} 0,{d_{im}} - C_{i'm}^d < {T_\delta }\\ {d_{im}} - C_{i'm}^d,{d_{im}} - C_{i'm}^d \ge {T_\delta } \end{array} \right. $

(6)

式中：δ_im表示特征点k_i的描述子残差δ_i中第m维的数值；i'表示特征描述子d_i根据码书编码的索引；T_δ是一个残差保留的判决门限，只有大于T_δ的残差才被保留，T_δ的引入使编码残差所用的比特数更少。最后，描述子编码码流、位置编码码流和低频子带编码码流共同组成了编码端输出的编码码流。

2 解码实现

解码过程的流程图如图 4所示。在解码端首先将低频子带进行JPEG2000解码，得到解码后的低频子带L'。然后进行高频子带的解码和复原。

2.1 高频子带解码

高频子带的编码信息包括特征描述子的位置信息和特征描述索引。对于位置信息，通过二进制游程解码进行无损还原。而描述子信息则通过索引找到码书中对应的描述子，并加上解码后的残差向量组成解码后的描述子。描述子的解码过程可表示为

$ {{d'}_i} = C_{i'}^D + {{\delta '}_i} $

(7)

式中：d_i'表示解码后的特征点k_i的描述子，δ_i'为解码后的残差向量。

然后将解码后的所有特征描述子逐一与库中描述子根据式(4)进行匹配。对于描述子d_i'，库中使dist_ij最小的描述子d_j^D对应的特征点k_j^D即为图像库中与特征点k_j^D匹配的特征点。匹配点k_j^D可以表示为

$ k_j^D = {\rm{M}}\left( {{k_i}} \right) $

(8)

为了使匹配更加有针对性，首先采用SPM (spatial pyramid matching，SPM)^[19]算法检索出与低频图像相似度最高的图像，形成图像库的一个子集，特征匹配只在图像库子集中完成。

通过这一过程，就可以建立一个从特征点集合K到库中特征点K^D的一个映射关系，K中所有元素都可以在K^D中找到一个对应元素M(k_i)。在图像库中，M(k_i)具有与k_i最为接近的特征描述子，说明M(k_i)特征点周围的图像块与k_i周围的图像块形态也十分接近。因此将k_i周围的图像块用于k_i周围图像块的复原。将从库中提取到的图像块进行透视变换，一个3×3的变换矩阵的基本形式为

$ \mathit{\boldsymbol{H}} = \left[ {\begin{array}{*{20}{c}} {{a_{11}}}&{{a_{12}}}&{{a_{13}}}\\ {{a_{21}}}&{{a_{22}}}&{{a_{23}}}\\ {{a_{31}}}&{{a_{32}}}&1 \end{array}} \right] $

(9)

变换矩阵H可以通过匹配特征点位置的对应关系求得。图像点的透视变换过程可以表示为

$ \left( {{x_i},{y_i},1} \right) = \left( {x_j^D,y_j^D,1} \right) \times \mathit{\boldsymbol{H}} $

(10)

将M(k_i)所在的高频图像进行n层小波分解，可以得到3×n个高频子带，M(k_i)在每个高频子带上都有一个对应位置，将对应位置周围的图像块经过透视变换后，分别贴补到待复原的对应高频子带中，贴补位置为k_i在待复原的高频子带上的对应位置。就可以得到复原后的高频子带。高频子带的复原过程如图 5所示。

2.2 小波逆变换

在完成低频子带解码和高频子带复原后，根据解码后的低频子带和复原后的高频子带，通过n层小波逆变换得到重构图像。

3 实验结果与数据分析 3.1 实验参数设置

本文的测试图像与图像库均来自于INRIA Holidays database图像库^[20]。这是一个自然图像库，库中共有1 490张图像，在这些图像中，有一些对相同场景从不同角度拍摄的图像。并另选取20张与库中图像场景相同，视角不同的图像作为测试图像。测试图像中，原始图像的分辨率分别为：1 600×1 200，4幅；2 048×1 536，8幅；2 560×1 920，6幅；2 816×2 112，2幅。

对于小波层数的选择，从提高压缩率的角度分析，小波变换的层数应该是越多越好，但是考虑到图像的复原效果，在实验中，当小波层数大于3层时，复原图像会出现比较明显的块效应，采用3层小波变换时，对于一些没有特征的区域，图像恢复效果不够理想。综合分析，本文采用2层小波变换。然后对小波变换后得到的低频子带进行JPEG 2000编码。

从高频子带图像中提取SIFT描述子，进行矢量量化编码，为了选取合适的码书大小，本文通过式(11)计算不同码书大小下，对待编码图像描述子预测的均方误差，结果如图 6所示。

$ \Delta \omega = \sum\limits_{i = 1}^n {\left( {\frac{{{{\left( {{d_i} - c_{i'}^D} \right)}^2}}}{{d_i^2}}} \right)/n - \omega \left( D \right)} $

(11)

式中：d_i和c_i'^D分别代表输入图像的描述子和与其对应的码书中的预测描述子，n为输入图像特征点的总数，ω(D)为用库中原始数据对输入图像描述子进行预测时的均方误差。

由图 6可知，随着码书的增大，描述子预测的均方误差会随之变小，当码书大小达到2¹³时，基本与未经矢量量化之前的预测结果接近，继续增大码书，预测均方误差趋于平稳，因此，本文选取的码书大小为2¹³。然后根据式(6)对描述子残差进行计算。

3.2 实验结果及质量评价

编码低频子带需要的码流大小与图像的复杂程度有很大关系，平坦区面积越大的图像，压缩低频子带所需的比特数越少，平均需要31.2 kB。根据图像的不同，每张高频图像中提取出的SIFT特征点个数差异很大，在测试图像中，特征点数最少100个，最多2 006个，平均833个。特征点数越多，编码特征点数所需的比特数就越多，编码特征的平均码流大小为3.61 kB。编码一个SIFT特征平均需要36 bit。计算每组实验的压缩比，可以得出，本文的编码方法平均压缩比可达437，最高可达1 388，这样的压缩比可以将一幅14 M的图像压缩到几十kByte甚至几kByte。

将压缩后的码流进行解码，并与JPEG 2000和JPEG进行对比。根据表 1中每张图像的压缩比，设置JPEG 2000与JPEG的压缩率，使本文算法的压缩比等于或略高于其他两种算法。然后在这个压缩比下，通过峰值倍噪比(peak signal to noise ratio, PSNR)与结构相似性(structural similarity index, SSIM)^[21]指标来对三种算法的性能加以衡量，PSNR是最经典的图像质量评价标准，而SSIM是一种新兴的图像评价标准，更符合人眼观看的视觉效果。

在表 1所示的压缩比下，分别比较本文算法与JPEG 2000、JPEG算法的复原图像效果，PSNR对比如图 7所示。

根据PSNR和SSIM的结果进行分析，将本文算法与JPEG算法进行对比，JPEG算法的压缩比在大于一定数值以后图像质量明显下降，有些图像无法达到与本文算法相同的压缩比，在此情况下，JPEG算法重构图像的PSNR值平均为26.63 dB，SSIM的平均值为0.774，本文算法所有重构图像的PSNR与SSIM数值均大于JPEG算法，本文算法的重构图像PSNR平均值为31.18 dB，SSIM的平均值为0.854。同等压缩比下，PSNR平均比JPEG高4.55 dB，SSIM平均比JPEG高0.08。用JPEG 2000算法对测试图像进行编码，编码一张图像需要的平均码流大小为35.42 kB，平均压缩比为418。本文算法的PSNR平均值比JPEG 2000高0.18 dB，SSIM平均值比JPEG 2000高0.001。以序号为6的图像为例，重构图像的主观效果的局部放大对比如图 8所示。为了更好地分析算法性能，在主观对比中，加入了相同压缩比下，与最新压缩标准HEVC帧内编码的对比。通过主观对比可以看出，JPEG算法的块效应较为明显，并且细节和色彩出现大量丢失。JPEG 2000算法对色彩重构较好，但一些细节部分会出现模糊，而本文算法对局部细节恢复较好。与HEVC算法相比，本文算法在颜色复原度上略优于HEVC。

4 结论

1) 本文提出了一种小波变换下基于SIFT特征匹配的图像编码方法。将空域图像转换到小波变换域，并在小波域完成特征匹配，压缩及图像重构。本文算法的平均压缩比可达到437。在相同压缩比下，重构图像的平均PSNR和SSIM均较JPEG与JPEG2000有所提升。

2) 由于编码过程加入了特征提取，本文算法的编码时间要长于JPEG和JPEG2000。但是解码过程不涉及到特征提取，因此本文算法可以实现快速解码。