由于天气、污染等原因导致的雾霾是一种普遍的现象。在雾霾天气，存在很多诸如PM2.5的大气粒子，它们不仅吸收和散射场景的反射光，而且还将一些大气光散射到相机。因此，相机接收到的图像严重降质，能见度大大下降^[1]。由于图像的质量下降，特征提取、目标跟踪以及目标识别等也受到干扰。因此，十分有必要研究图像去雾算法，提高视觉系统的环境适应性。

随着计算机技术、光学器件技术的提高，图像去雾算法受到了广泛的关注，取得了很多成果，也应用在许多领域。常见的去雾算法可以分为两类：单幅图像去雾和多幅图像去雾。Raanan Fattal^[2]优化了大气成像模型和大气衰减模型，加入了表面阴影量，提高了清晰度和对比度。Kaiming He^[3]提出了暗通道先验的算法，提升了算法的速度以及天空区域的识别率。此外还有很多学者^[4-6]基于传统算法提出了改进方法，但是基本的原理都与传统算法相同。随着神经网络的发展，Cai B^[7]基于深度学习提出一种名为DehazeNet的可训练的端到端系统，用于参数估计，去雾效果相较于传统算法，色彩还原度高，清晰度高，但是速度慢，需要前期大样本训练。

对于单幅图像去雾，使用这些方法可以有效地改善模糊图像的质量，但是由于需要许多未知参数，因此它们通常依赖于输入图像的某些假设或先验知识^[8]。但是多幅图像去雾算法不存在这样的问题，通过输入的多幅图像，可以轻松地获得去雾后的图像，具备较强的去雾能力。在多幅图像去雾算法中，偏振去雾算法是一种兼备速度快、效果好以及鲁棒性强等优点的算法^[8-12]。

本文针对现有偏振去雾算法鲁棒性不强和图像增强效果有限的问题，提出了一种基于多尺度奇异值分解的偏振图像融合去雾算法，并设计了不同能见度条件下的验证实验。与传统偏振去雾算法相比，本文算法具有较强的自适应性和鲁棒性，能够有效地改善传统算法中出现的光晕以及天空区域过曝的问题，提高了复杂光照条件下的视觉成像质量，具有较大的应用前景。

在雾霾天气下，大气中悬浮大量的微小颗粒，大气光和目标反射光会与粒子发生散射。雾霾条件下相机成像模型如图1所示，成像系统接收光由两部分构成：由大气中的粒子对于光源(${{{A}}_\infty }$)的散射之后到达成像系统的大气光(${{A}}$)；目标反射光强(${{L}}$)经过大气中的粒子吸收和散射后的直接透射光(${{D}}$)。

图  1  雾霾天气下的相机成像模型

Figure 1.  Block diagram of the image measuring system

则未经大气散射衰减和没有大气光干扰的目标反射光强为^[13]：

上式中总光强(${{I}}$)可以由相机直接获得，通过估算无穷远处的大气光(${{{A}}_\infty }$)以及大气光(${{A}}$)，即可计算出理想的目标反射光强，达到图像去雾的目的。

传统的偏振去雾算法中^[8-12]，通过旋转偏振片获得最亮和最暗的偏振图像，利用得到的图像进行差分求取(1)中的参数。然而，这种方法需要人为主观的判断最好与最坏的图像，导致了算法鲁棒性不强，不易操作，并且人为主观的判断也存在较大的误差。因此，Jian Liang^[12]等人提出了基于Stokes矢量的去雾算法，首先在相机上安装偏振片，分别获得0°、45°、90°和135°的图像，Stokes矢量可以表示为^[14]：

在（2）中${{{S}}_{\rm{0}}}$表示入射光的总强度，${{{S}}_{\rm{1}}}$和${{{S}}_2}$表示偏振态的入射光。由于影响图片质量的为白噪声，四幅图像相加以及两幅图像的相减均可以达到降噪的效果。在公式(3)中计算了偏振角(${{\theta }}$):

根据文献[15]，大气光可以表示如下：

其中

式中${{{\theta }}_A}$为出现在整幅图像中出现概率最大的偏振角。

无穷远处的大气光(${{{A}}_\infty }$)可以通过暗通道先验的方法获取^[3]。暗通道图像可由下式获得：

选取单个像素在红黄蓝三个通道的最小值作为该像素的亮度值，然后选择合适的窗口大小${\bf{\Omega }}(x)$，在该窗口中取所有像素的最小值作为窗口所有像素的亮度值，即可得到原始图像的暗通道图像。在[3]中，作者根据实验经验选取窗口大小为15×15。在[16]中作者实验发现对近处物体取大窗口，远处物体取小窗口重构效果好的特点，并依据估计的深度设计了不同的窗口大小。在本文中，根据实验结果，选取窗口为15×15进行计算。

得到暗通道图像后，选取亮度超过设定阈值的像素求平均值，得到无穷远处的大气光(${{{A}}_\infty }$)。将（4）中求取的大气光(${{A}}$)和无穷远处的大气光(${{{A}}_\infty }$)代入（1）式即可求得去雾后图像。

经典的Schechner^[9]算法以及2.1节使用的方法均假设大气光是部分偏振光，直接透射光是非偏振光。但是目标反射光的偏振特性不仅与大气中的微粒有关，还与目标本身的材料、粗糙度以及颜色有关。Fang Shuai^[15]拍摄了同一场景的5000张图像，证明了目标反射光的偏振特性不可以忽略，并提出了改进的雾天偏振成像模型：

在上式中，总光强${{I}}$和总偏振度${{P}}$可以直接通过原图像计算得出。所以在该模型中，只需要估计出目标光偏振度${{{P}}_D}$，大气光偏振度${{{P}}_A}$以及无穷远处的大气光强${{{A}}_\infty }$即可实现去雾。

合肥工业大学的高隽^[16]通过高斯滤波求取大气光。一般的大气光(${{A}}$)估算方法多是针对非偏振图像进行滤波处理，这类方法没有考虑深度对于偏振度的影响，导致距离视点较远区域重构效果较差。因此，高隽提出使用高斯滤波对大气光偏振度${{{P}}_A}$进行估计。

高斯滤波是一种根据正态分布来分配权值的线性平滑滤波器。相对于均值滤波器和中值滤波器，高斯滤波器引入了空间距离因子，对滤波窗口内的像素计算权重，越位于窗口外侧的像素权重越低，其对待处理像素的影响就越小，位于窗口内侧的像素与待处理像素更接近，权重更高。高斯滤波数学表达式如下：

通过相机获取四个偏振角度（0°、45°、90°、135°）下的偏振图像进行（9）中的高斯滤波处理，得到不同偏振角度下的大气光估计图像$GB({{{I}}_0})$、$GB({{{I}}_{45}})$、$GB({{{I}}_{90}})$、$GB({{{I}}_{135}})$。

式中，${{\omega }}$为权重矩阵，根据去雾的效果需要人工调节，x和y表示周边像素对应于中心像素坐标的相对坐标，$\vartheta $表示模糊半径即权重矩阵的大小；${{I}}$表示待滤波图像；${\bf{\Omega }}$为窗口大小，对于不同的图像，具有不同的最优窗口，可以根据去雾效果好坏，依据经验或者手动调节选取窗口大小。

使用滤波得到的大气光估计图像计算大气光偏振度${{{P}}_A}$，由${{{P}}_A}$可求得直接透射光强(${{D}}$)以及目标光偏振度${{{P}}_D}$^[15]。将求得的参数代入模型(7)中，即可得到去雾后的图像${{L}}$。

立足于提高传统算法的适应性，提高去雾图像的质量，本文设计了如图2所示的去雾算法流程。首先，使用了基于最小二乘的方法计算了更加精确的偏振信息，改善了以往偏振信息计算不准确的问题。然后，分别使用基于高斯滤波的去雾方法以及基于Stokes矢量的去雾方法处理图像。最后，根据两种算法去雾结果的互补特性，选择计算量较小、对于高像素图像处理更简便的图像融合算法，对两种去雾结果进行融合，得到去雾效果更好的图像。

图  2  多尺度奇异值分解的偏振图像融合流程图

Figure 2.  Flowchart of multi-scale singular value decomposition polarization image fusion

现有的偏振去雾算法本质上都是利用Stokes矢量的各分量直接求偏振角和偏振度，这种方法存在求取偏振信息噪声大，区分度不强等问题。目前可以通过偏振相机直接获得四个方向的偏振图像，使用各个分量分别计算偏振信息实质上并没有充分利用冗余的信息。因此，本文使用最小二乘的方法计算偏振信息，经过实验验证，基于最小二乘计算偏振信息较传统方法噪声小，物体辨认更加清晰。

如图3是基于Sony IMX250芯片的Blackfly偏振相机，相机可以通过单次曝光，自动实时地对全天域偏振模式进行测量，从而获取天空的偏振信息。

图  3  Blackfly型偏振相机

Figure 3.  Blackfly Polarization Camera

传感器具有若干偏振测量单元，如图3(b)所示，每个测量单元由4个像素组成，对应像素具有不同偏振敏感方向，分别为0°，45°，90°和135°，并按重复的两个像素块排列。每个像素的偏光滤镜（C）均涂有抗反射层（B），并位于镜头（A）和感光光电二极管（E）之间。

在每个像素中，忽略足够小的圆偏振光，则CCD相机主要测量线偏振光。相机的任一像素的光强响应可以表示为^[14]：

其中，${{{I}}_x}$为像素的亮度值，${{{K}}_x}$和${{{\alpha }}_x}$分别表示相机感光系数和雕刻光栅的极化角度，可以通过标定获得。

所求未知数共有三个${{I}}$、${{P}}$和${{\theta }}$，而Blackfly型偏振相机可以方便地获得四个方向的偏振图像，则重新展开(10)式可得:

利用偏振相机得到的偏振图像联立方程可得：

其中，${{{\rm A}}} = \left[ \begin{gathered} \cos 2{\alpha _1}{\rm{ }}\sin 2{\alpha _1}{\rm{ 1}} \\ \cos 2{\alpha _2}{\rm{ }}\sin 2{\alpha _2}{\rm{ 1}} \\ \cos 2{\alpha _3}{\rm{ }}\sin 2{\alpha _3}{\rm{ 1}} \\ \cos 2{\alpha _4}{\rm{ }}\sin 2{\alpha _4}{\rm{ 1}} \\ \end{gathered} \right]$，${{Q}} = \displaystyle\frac{1}{2}\left[ \begin{gathered} {{IP}}\cos 2{{\theta }} \\ {{IP}}\cos 2{{\theta }} \\ {\rm{ }}{{I}} \\ \end{gathered} \right]$， ${{U}} = \left[ \begin{gathered} {{{I}}_1}/{{{K}}_1} \\ {{{I}}_2}/{{{K}}_2} \\ {{{I}}_3}/{{{K}}_3} \\ {{{I}}_4}/{{{K}}_4} \\ \end{gathered} \right]$。

最小二乘估计为：

则偏振角$\theta $和偏振度$P$可由下式计算：

其中，${{{q}}_1}$、${{{q}}_2}$和${{{q}}_3}$为${{Q}}$的向量元素。

图像融合可以分为像素级、特征级以及决策级融合^[17]。多分辨率奇异值分解和基于小波变换的图像融合是像素级融合方法的重要构成部分。对比基于小波变换的图像融合，基于多分辨率奇异值分解的图像融合性能与之接近，但是计算更简单、实时性更加突出以及对复杂、高像素图像处理简单方便。因此，在本文中使用基于多分辨率奇异值分解的图像融合。

多分辨率奇异值分解源于Burt和Adelson提出的拉普拉斯金字塔算法^[18]，类似于小波变换，基本思想是在平滑分量的每一层上用奇异值分解代替滤波。图4展示了三层MSVD分解的结构。假设${{X}}$是一个$M \times N$矩阵，将${{X}}$分成不重叠的${\rm{2}} \times {\rm{2}}$块，并且通过堆叠，以形成矩阵${{{X}}_1}$，并将每个块排列成${\rm{4}} \times {\rm{1}}$的向量。${\rm{4}} \times {\rm{4}}$散布矩阵的特征分解是：

将（16）中所得的奇异值以降序来排列：

令${{{\hat X}}_1} = {{U}}_1^T{{{\bar X}}_1}$，${{{\hat X}}_1}$第一行包含对应于最大特征值的近似或平滑分量，剩余行包含对应于最小特征值的细节分量。每行中的元素可以被重新排列形成一个$M/2 \times N/2$的矩阵，表示为${{\bf{\Phi }}_1}$。以类似的方式，将剩余三行${{{\hat X}}_1}(2, \cdot )$、${{{\hat X}}_1}(3, \cdot )$、${{{\hat X}}_1}(4, \cdot )$，重新排列成M/2×N/2矩阵，分别表示为${\bf{\Psi }}_1^V$、${\bf{\Psi }}_1^H$和${\bf{\Psi }}_1^D$。下一层变换将${{\bf{\Phi }}_1}$替换成${{X}}$，重复上述步骤。整个分解过程可以描述如下：

图  4  三层MSVD分解结构示意图

Figure 4.  Diagram of the three-layer MSVD decomposition structure

被融合的偏振去雾图像${{{I}}_1}$和${{{I}}_2}$分别用MSVD分解为$l$(${\rm{l}}=1,2,\cdots,L$)层。在每个分解层次中，细节分量系数对应于图像中更尖锐的亮度变化，如边缘图像等。而较大的奇异值能更好地保持原始低分辨率图像的有用信息，因此融合规则将选择两个MSVD细节分量的较大值。对于底层($l$=1)，融合规则取MSVD平滑分量系数的平均值，因为底层的平滑分量系数是原始图像的平滑和亚采样版本[17]。融合图像可通过下式得到:

实验整体设备和拍摄场景如图5所示，偏振相机通过USB3.0接口与计算机相连，使用支架固定偏振相机保持其稳定性。

图  5  实验装置图

Figure 5.  Experimental setup

调节相机的焦距以及光圈，使得相机尽可能成像清晰。通过一次曝光，即可得到偏振角度分别为0°、45°、90°以及135°的偏振图像。在不同环境条件采集图像时，温度、气压、湿度以及能见度等均对成像有一定影响，但是能见度才是决定相机成像好坏的主要因素。在本文实验设计时，主要对比了几种不同能见度条件下算法的去雾效果，以验证算法的适应性。

为检验改进的偏振信息计算方法效果，使用相机拍摄不同条件下的偏振图像。如图6(a)为晴朗无雾霾天气拍摄的图像，拍摄时能见度30千米；6(b)为雾霾天气拍摄图像，能见度为4千米。

图  6  偏振信息计算结果示意图

Figure 6.  Calculation results of polarization information

由于目标的偏振特性与目标的本身的材质、颜色以及光滑度有关，所以天空与建筑物、不同建筑物之间以及建筑物与树木之间的偏振信息肯定是存在差异的。只要偏振信息的计算方法不存在错误，那么偏振信息的计算结果越清晰，噪声越小，则说明计算的结果越精确，误差越小。

使用3.1节中方法计算其偏振角图像，分别得到偏振角和偏振度计算结果。从图6(c)可以看出，在晴天的偏振角图像中，建筑物与天空的区分度不强，计算得出的建筑物偏振角与天空偏振角差异很小，难以区分天空以及建筑物。在图像下部分建筑物的轮廓模糊不清，几乎被噪声覆盖。整幅图像计算的偏振角度差异不明显，噪声很大。在雾天的计算的结果中，除了边界部分，建筑物与天空的偏振角也基本相等。与晴天的计算结果类似，建筑物的轮廓被噪声掩盖，难以区分。因此，传统的计算方法计算偏振角信息时引入了大量的噪声信息，存在缺陷。

如图6(e)是使用最小二乘方法计算的晴天偏振角度，6(f)是使用最小二乘方法计算的雾天偏振角度。无论是晴天还是雾天，相比于传统方法，可以看出建筑物的轮廓更加清晰，与天空的区分度也更好，引入的噪声更小。使用信息熵，平均梯度以及原图与结果的相似性(SSIM)评价偏振信息计算结果：

由表1可知，在晴朗条件下，最小二乘法的三项指标分别比传统方法高出14%、17%以及0.1%。因为SSIM指标已经与1非常接近，所以提高幅度较小，平均梯度以及信息熵两项指标均有明显的提升。

由表2可知，在雾天条件下，最小二乘法的前两项指标分别比传统方法高31%、15%。从以上结果来看，无论是晴天还是雾天，最小二乘方法的绝大多数性能指标均优于传统的计算方法。因此采用最小二乘代替传统的偏振信息计算方法，可以提高偏振信息的计算精度，降低噪声。偏振信息计算的越准确，最终的去雾的效果也将更好。

为了验证图像融合算法的可行性，使用相机采集了不同能见度条件的偏振图像，使用3.2节所述的算法流程处理图像。

如图7（a），拍摄时能见度为3.9千米，从图中可以看到，远处的建筑物被雾霾遮挡，近处的建筑物也模糊不清，能见度大大下降，图像的细节丢失。

图  7  图像融合去雾流程示意图

Figure 7.  Schematic diagram of the image fusion defogging process

如图7（b），是使用基于高斯滤波的去雾算法结果。在实际操作中，使用原作者提供的参数，并没有进行过多的调节。去雾后整体图像的色调偏暗，图像有一些失真，建筑物的细节丢失，参数的适应性较差。但是，通过观察可以发现图像去雾有一定的效果，原本远处模糊的建筑物轮廓有所增强。

如图7（c），是使用基于Stokes矢量算法的去雾结果，相比于原图，总体上去雾的效果较好，远处的建筑物能够清晰辨认，近处建筑物的细节也有所加强。但是不难发现，该图像在天空区域以及图像下方的白色雷达区域出现了过曝的现象，并且在建筑物的边缘周围还出现了光晕现象。过曝的问题通过调节模型中的参数，可以得到一定程度改善，但是实验发现，很难找到最优的参数，并且往往过曝现象只能得到改善而不能消除。

从以上分析可知两种方法的优缺点存在互补，因此可以采用图像融合的方法对两种结果的优点加以结合。在图7（d）中展示了两种算法结果的融合结果。从图中可以发现，与基于Stokes矢量的方法相比，天空区域的过度曝光现象已经消除，改善了建筑物的光晕现象与白色雷达的过度曝光现象。与基于高斯滤波的方法相比，丰富了图像的细节，图像更加接近平时所观察到的实际场景。此外，相比于两种算法，图像融合可以代替人工调节参数，增强了算法的鲁棒性。

为了验证算法能够在不同的能见度条件下的适应性，分别在能见度为6.8千米（图8）和9千米（图9）时做采集数据，并按照所提出的算法进行数据处理。编号为a的为采集的原图，编号为b的为基于高斯滤波的去雾结果，编号为c的为基于Stokes矢量的去雾结果，编号为d的为图像融合的结果。

图  8  能见度6.8千米去雾结果

Figure 8.  Defogging results with a visibility of 6.8 km

图  9  能见度9千米去雾结果

Figure 9.  Defogging results with a visibility of 9 km

通过观察去雾结果可以发现，在能见度较高的情况下，相比于现有的两种算法，使用的算法能在一定程度上提高去雾效果，增强了细节、改善了过曝现象。在能见度较低的情况下，提出的算法可以有效的改善过曝现象，显著的提高图像质量，丰富了图像的细节信息，达到了去雾的效果。

为对比算法的去雾结果，本文使用开源基于MSCNN^[19]的去雾算法分别对三种可见度的图片进行了处理，结果如图10所示。

图  10  MSCNN去雾结果

Figure 10.  Defogging results of MSCNN

从图10结果来看，基于MSCNN的去雾算法有一定的去雾效果，但是同时在天空区域出现了不同程度的过曝现象。在表3中，对比了几种去雾算法的信息熵以及平均梯度。

此外，本文算法在i7-8500U 1.80Ghz计算机上使用Matlab2017b上实现，表4中记录了四种去雾算法处理一副图像的计算时间。

在不同的天气条件下，基于Stokes矢量的方法信息熵一般较大，因为基于Stokes矢量的方法去雾效果较好，细节较为丰富，但是出现了光晕以及过曝的现象，影响了图像的整体质量，降低了图像的平均梯度。基于MSCNN的去雾算法有一定的去雾效果，提高了图像的信息熵以及平均梯度，但是在天空区域也出现不同程度的过曝现象。基于高斯滤波的方法在信息熵以及平均梯度上都不是最优的，但是图像没有出现光晕与过曝现象。

由表3可知，图像融合提高了雾霾天气的信息熵，在三种天气条件下，对比于原图，信息熵最大提高了1.3%，相比于高斯滤波最大提高了18.9%。此外，无论在什么天气条件下，图像的平均梯度通过图像融合都得到了显著的提升，相比于基于高斯滤波算法、基于Stokes矢量算法、以及基于MSCNN算法最大分别提高了27.4%、38.4%和20.4%。

从算法的时间成本来看，本文的算法相较于传统去雾算法增加了图像融合的计算步骤，但是基于多分辨率奇异值分解的图像融合具有计算简单，对复杂、高像素图像处理更方便的优点。根据表4各算法时间成本的对比可知，图像融合步骤的时间消耗仅为基于Stokes矢量算法的23.3%，相比于基于滤波的算法来说，也没有显著增加时间消耗，总耗时与基于MSCNN的去雾算法持平。因此，本文提出的算法在丰富了图像的细节，改善了光晕现象的同时，没有大幅度增加整个去雾流程的计算时间。

本文针对现有偏振去雾算法鲁棒性不强和图像增强效果有限的问题，提出了一种基于多尺度奇异值分解的偏振图像融合去雾算法。采用最小二乘法，提高了偏振信息测量的准确度；在此基础上，给出了一种基于多尺度奇异值分解的偏振图像融合去雾算法，并设计实施了不同能见度条件下的验证实验，结果表明，在没有大幅增加时间成本的条件下，图像信息熵与平均梯度最大可分别提高18.9%和38.4%。采用图像融合的方法具有较强的自适应性和鲁棒性，能够有效地改善经典去雾算法中出现的光晕以及天空区域过曝的问题，对于利用可见光的偏振特性进行图像去雾和图像增强效果较好，具有较强的工程实用性。总体来说，该算法有效地提高了复杂光照条件下的视觉成像质量，在交通路况监控、航拍、军事目标探测与识别等视觉系统的应用及发展前景广阔，具有重要的研究意义。