融合Retinex和离散小波奇异值分解的远距离目标图像清晰化

徐兴贵，杨润华，冉 兵，樊香所

（1.云南财经大学 信息学院，云南 昆明 650221；2.中科院光电所 自适应光学重点实验室，四川 成都 610209；3.宜宾学院 智能制造学部，四川 宜宾 644000）

引言

有别于常规成像场景，大气湍流和气溶胶等光传输介质[1]对光的吸收和散射影响严重，使得远距离成像系统获得的目标图像同时出现噪声、模糊和低对比度等退化降质。因此，针对卫星遥感成像、天文光学系统成像和远距离可见光或红外成像系统等场景下的有效低照度降质图像清晰化（增强）面临较大挑战。

目前，图像清晰化方法大致可分为7 类：直方图均衡法[2-4]、基于Retinex 的增强方法[5-7]、基于物理模型的复原方法[8-9]、基于融合模型的复原方法[10-11]、深度学习的清晰化方法[12-13]、多应用背景[14-17]及图像处理下[18-20]基于直方图的规定化增强方法[21]和基于灰度值的拉伸法[22]。直方图均衡法利用图像像素值的直方图分布特征指导增强图像，如广义自适应直方图均衡[2]、奇异值均衡[3]以及照度直方图均衡[4]。由于均衡策略依据直方图统计特性，此类均衡方法容易造成过增强或者欠增强，并且不会有去模糊的效果。图像拉伸主要用来改善图像对比度，图像拉伸可具体分为3 种方式：基于图像像素的灰度拉伸增强法，基于图像像素邻域关系构建的直方图均衡化增强方法和基于变换的直方图规定化增强。特别地，直方图规定化方法能够有选择地控制直方图变成规定的形状，从而弥补了一般均衡法不具备的交互特性。第2 种基于Retinex 的图像增强方法则基于颜色恒定性假设，即物体表面的反射特性决定图像质量，与成像系统中的光辐射分量不相关。该类方法特别适合光照度不充分的成像场景。近年来，Rahman 等人[6-7]提出多尺度Retinex 算法提取物体反射分量，进而提取出包含图像细节信息的清晰本征值获得增强效果。另一种基于物理模型的清晰化方法则集中在物理模型的构建和先验信息的利用上。例如：Tarel 等人[8]提出的大气散射模型法，He 等人[9]提出的暗原色先验理论找出图像局部暗原色信息并据此优化透射率后得到清晰图。此类方法利用先验信息在常规自然成像场景中取得了较好效果，但远距离成像场景显然不满足条件。此外，融合模型的方法则采取序列“串联式”的策略增强图像，如Ancuti 等人[11]提出的图像融合方法，首先从退化图像中估计退化权重谱，接着根据权重谱推导介质散射和吸收模型，最后实现图像清晰化。然而，该方法的权重谱估计依赖大量工程经验，泛化能力差。融合方法中，小波变换已被广泛应用于图像处理，如特征提取[16]、去噪[17]、图像压缩和超分辨率[18-19]等，取得了许多较好的效果。最近，基于深度学习的方法被广泛用于图像清晰化，例如，Chen 等人[12]则利用R-CNN网络模型来处理CT 图像，Yao 等人[13]提出一种基于概率模型的深度学习方法清晰化降质图像。然而，该类方法的深度学习网络优化和数据标注等方面的软硬件成本依然较高。值得一提的是，远距离降质图像的高清参考图像以及光传输物理成像模型先验信息等难以获取，导致基于物理模型和深度学习的方法难以实现。据此，本文采用Retinex和小波变换的融合模型实现远距离低照度降质图像增强。

结合远距离成像场景下的低照度和多种退化降质图像特点，受Retinex 恒常性原理和小波重构高分辨图像的启发，本文采用融合离散小波低频分量的奇异值分解策略来清晰化远距离退化目标图像，并对其有效性进行了验证。所提方法主要分为两步：首先利用Retinex 方法“粗”提取降质图像的照度分量和反射分量，然后利用小波奇异值分解“精”提取反射分量，并通过小波逆变换实现图像清晰化。

1 所提模型

1.1 照度分量提取

如引言所述，Retinex[5-7]增强模型基于颜色恒定性假设，能够将退化图像分解为低频照度分量和变化剧烈的高频反射分量。Retinex 模型表达为

式中：Ii为退化降质图像；λ为信道波长；Li和Ri分别为照射分量和反射分量。估计光照分量和反射分量则可通过环绕尺度函数滤波法[15]，照度分量表示为

以及

式中：Fi(x,y)表示环绕尺度函数；⊗为卷积操作算子。相应地,反射分量R（x,y）可以通过（1）式的对数变换解算。一般地，滤波器尺度参数 σi的值越大，照度分量估计的全局性越好，但局部性越差。为兼顾全局性和局部性，本文采用作者曾提出的自适应全尺度AFSR 方法[15]估计照度分量和尺度分量，其尺度参数映射公式为

图1 AFSR 分解结果Fig.1 Results of AFSR decomposition

1.2 小波子带奇异矩阵分解

二维小波变换将图像分解成不同的频率范围，与文献[14,16]中方法相似，本文采用离散小波变换（discrete wavelet transform,DWT）先沿图像的行进行一维离散小波变换，然后进行图像列分解，由此产生4 个类型的子带图像，即称为低-低（LL）、低-高（LH）、高-低（HL）和高-高（HH）。值得一提的是，LL 小波变换子带对应为图像低频照度分量，进一步通过奇异值分解（singular value decomposition,SVD）提取出退化图像的高频分量，最后经过反向离散小波变换（inverse discrete wavelet transform,IDWT）得到清晰图像。奇异值分解可实现图像信号的稀疏分解，SVD 矩阵分解可由前r个主奇异值和特征向量近似构成原始矩阵A，表达式为

式中T 为转置操作符。原始图像矩阵可分解为特征矩阵U、奇异值构成的对角矩阵Σ和系数矩阵V。一般地，奇异值σi按降序排列构成对角阵，其前面10%甚至1%的奇异值之和就占了全部奇异值总和的99%。通常，采用生成的归一化矩阵（均值为0，方差为1）最大奇异值便可求取主奇异值，修正系数公式可表达为

结合（5）式和（6）式，原始矩阵经SVD 分解后可表达为

进一步，新的LL 分量的奇异值分解表达为

通过逆离散小波变换，最终的清晰化图像Aobj可通过分解后的LL 子带和其他3 个子带重建，即：

实验中，我们使用了db.9/7 小波函数作为小波变换的母函数。如图2所示，离散小波变换分解原退化图像得到4 个子波带HH/HL/LH/LL。可以看出，低频子波带图LL 刻画了原图像灰度缓慢地变化的大部分信息，代表了图像的大致概貌和轮廓，是原图像（图7 星云图）的近似信息。相应地，低-高（LH）、高-低（HL）和高-高（HH）3 个子波带刻画的是灰度变化较快的区域，表现为图像边缘区域的部分。从图2 中可以看出，HH 子波带的边缘区域部分最为明显。

图2 离散小波分解所得的4 个子波带Fig.2 Four wavelet subbands obtained from discrete wavelet decomposition

另一方面，如图3（a）和图3（b）所示，AFSR 增强后的LL 子波带奇异值主要集中在值较大的前半部分，而原退化图的低频LL 奇异值矩阵则比较均匀地分布，其原因是低频的照度分量模糊了图像，进而导致光电探测器靶面上的能量即灰度值均匀分布在各区域。据此，低频子波带的SVD 分解后，通过去除低奇异值对应的模糊分量信息以提高图像对比度。进一步，针对小波变换后获得的含有模糊信息的低频子波带，利用SVD 分解的方式去除冗余的模糊信息。如图3（c）和图3（d）所示，LL 子波带经过SVD 分解后，图像边缘轮廓区域得以保留而平滑区域的模糊信息被去除，最后通过小波逆变换后便可实现图像清晰化。

图3 奇异值矩阵图示和LL 子波带SVD 分解前后对比图Fig.3 Schematic diagram of singular value matrix and comparison of LL wavelet subbands before and after SVD decomposition

对比小波变换处理前后的降质卫星图像（图6(a)）的直方图和方差值如图4所示。可以看出,所提的融合算法能够大大提高图像的方差，即客观反映对比度提高，增强后的图像直方图得到了改善，更多评价指标见第2.1 节。

图4 所提融合算法处理前后结果对比Fig.4 Results comparison before and after processing of proposed fusion algorithm

1.3 算法流程

综上所述，可得出本算法结构框图。如图5所示，算法步骤归纳为

图5 所提算法结构框图Fig.5 Block diagram of proposed algorithm structure

步骤1）AFSR 分解得到反射分量和照度分量；

步骤2）离散小波变换分解得到4 个子波带HH/HL/LH/LL；

步骤3）根据（8）式计算SVD 系数 ξLLA¯；

步骤4）根据（10）式得到增强的图像。

2 实验结果与分析

本文实验平台为64 位win10 操作系统，Intel Celeron（R）CPU1.83 GHz、8G 内存，MATLAB2014b，图像资源来自文献[15,19,20]，原始图像如图6（a）、6（b）和6（c）所示，3 幅图分别代表典型的空对地卫星遥感图像、近地面远距离红外图像目标和地对空可见光星云图像，该3 种远距离成像条件下的图像均出现不同程度的模糊、噪声和低对比度降质。3 幅原始图像大小分别为375×327 像素、330×329 像素、335×330 像素。

图6 原退化降质图Fig.6 Original degraded images

本文利用极大化的方式获取经过离散小波变换后的低频分量LL 的奇异值修正系数，即所提方法中的低频子带系数采用极大值的选择方式。同时，增加了极大值的选择和基于区域方差选择两种方式的清晰化对比结果。如图7所示，利用图像处理后的对比度（纵坐标）刻画清晰化效果，5 种远距离目标图像场景下获取的结果表明所选择的极大值方式能得到更好的奇异值修正系数。

图7 低频子带系数修正对比Fig.7 Comparison of low frequency subbands coefficient correction

2.1 评价指标

远距离退化降质图像的清晰参照图像不易获取，故采用无参考的图像质量评价[1,15]指标，主要包括均熵（mean entropy,ME）、标准差（standard deviation,SD）、边缘强度（edge intensity,EI）和图像对比度指标（image contrast,IC）。均熵为无量纲的信号随机统计量，该值越大则表示图像信息量越大，也即图像边缘纹理等信息就越丰富。均熵ME 表达式为

图像标准偏差反映了图像的高频分量，该值越高，对比度越高，图像细节信息也越丰富。图像标准方差SD 表示为

边缘强度采用Sobel 算子来刻画图像边缘细节程度，边缘强度EI 可表示为

式中Sobel 算子 ∇ij根据文献[19]计算得到。

图像对比度指标（image-contrast metric,IconM）则将图像分割成16×16 块，得到每个块的最大和最小像素值，计算每块中的相对对比度并求和，该指标通过（14）式计算：

2.2 实验结果

图9 文献[3]算法的增强结果Fig.9 Enhancement results of algorithm in reference [3]

图10 文献[4]算法的增强结果Fig.10 Enhancement results of algorithm in reference [4]

图11 文献[19]算法的增强结果Fig.11 Enhancement results of algorithm in reference [19]

实验结果中，图8～图12分别对应文献[2]、文献[3]、文献[4]、文献[19]和所提算法的清晰化结果,表1、表2 和表3 分别为图6（a）、6（b）、6（c）3 个图在相应算法中的客观评价指标。从主观视觉度和清晰度看，文献[2]的算法仅次于所提算法，但图像的高亮度（偏白色）区域出现了过增强，且目标边缘位置出现了大量不光滑的毛刺，如图8（a）中上半部分的操场区域全部被增强为纯白色的区域，图8 中的增强图像出现了大量的噪声；文献[3]和文献[4]两种方法的增强效果最不明显，增强后的图像对比度不高，去噪效果不明显，模糊现象严重；文献[19]的增强效果一般。相比之下，所提方法在视觉清晰度和对比度方面均得到了较好的提高，如图12（c）星云图的纹理得到了复原，图12（a）右下角的房子和围栏等目标突出的部分也得到较好增强，图12（b）的红外模糊目标图像去除模糊效果和噪声最好，视觉清晰度最高。

图8 文献[2]算法的增强结果Fig.8 Enhancement results of algorithm in reference [2]

图12 本文算法的增强结果Fig.12 Enhancement results of our proposed algorithm

表1 图6（a）的清晰化结果评价Table 1 Enhancement evaluation results of Fig.6(a)

表2 图6（b）的清晰化结果评价Table 2 Enhancement evaluation results of Fig.6(b)

表3 图6（c）的清晰化结果评价Table 3 Enhancement evaluation results of Fig.6(c)

从客观评价结果看，本文算法的清晰化结果在图像信息熵、平均梯度、对比度和边缘密度等方面总体上都明显优于文献[3]和文献[4]的算法，值得注意的是，文献[3]的对比度值虽然比所提算法高，但其出现了不是我们所希望的欠增强结果。另外，文献[2]算法处理图6（b）的边缘密度指标最高，但究其原因主要是该算法结果中依然存在大量噪声。总体来看，所提算法在客观指标的综合评价方面最好。

3 结论

本文针对远距离成像系统获取的低照度降质图像清晰化问题，提出了一种融合Retinex和离散小波奇异值分解的图像对比度增强算法。该方法能够利用自适应全尺度Retinex（AFSR）方法“粗”提取照度分量和反射分量；同时，通过离散小波变换能够估计出低频子带图像奇异值矩阵并应用逆小波变换实现“精”重建图像。对比其他方法结果表明，所提方法在主客观评价指标方面均优于其他算法。下一步考虑将所提图像增强后处理方法扩展应用到远距离暗弱目标检测和识别。