开放科学（资源服务）标志码（OSID）：

0 概述

在恶劣的天气条件下获取户外图像时，捕获的场景能见度通常较低，这是由于光与悬浮在大气中的粒子发生相互作用，如散射、吸收、反射等，导致图像对比度和饱和度降低，颜色褪色，这会对许多户外计算机视觉应用（如识别、检测、跟踪、分类等）产生不利影响^[1-2]。因此，对图像去雾进行研究在图像处理和计算机视觉领域具有重要的意义和价值，设计有效的去雾方法是近年来的研究热点。

早期图像去雾研究采用图像增强技术来提高模糊图像的可见度并改善色彩。其中，Retinex^[3]和choi^[4]2种图像处理技术是典型代表。但是，基于图像增强的去雾方法由于没有考虑图像的退化机理，不能完全去雾。因此，后续研究工作重点关注基于大气散射模型的去雾方法，该类方法首先构建雾天成像模型，然后估计未知参数，从图像退化的原因出发对物理模型进行反演，从而得到无雾图像，但是，这些方法大都需要利用先验知识或假设来估计关键参数。例如，TAREL等^[5]使用中值滤波器来估计耗散函数，由于中值滤波器显示的边缘保留性能较差，因此该方法在图像景深变化处有雾残留。NISHINO等^[6]采用贝叶斯后验概率模型充分挖掘图像中潜在的统计特征从而进行去雾处理，该方法在浓雾图像上具有一定的优越性，但对于薄雾图像，其去雾后的效果颜色过亮，缺乏真实感。HE等^[7]基于对室外无雾图像的统计，提出暗原色先验去雾方法，该方法去雾效果明显，但对于天空区域，其易出现失真现象。GIBSON等^[8]在HE等^[7]研究成果的基础上提出中值暗通道先验方法，该方法不需要对透射率进行精细化，在一定程度上加快了去雾过程，但是，这种方法并不能获得良好的视觉效果。文献[9-11]检测天空区域并根据其雾度对检测区域中的透射图进行补偿，但是，由于很难确定雾度和天空区域，因此其结果并不准确。ZHU等^[12]建立一个线性模型来模拟色彩衰减下的场景深度，并用监督学习方法学习模型的参数，但是，由于大气散射模型中的散射系数实际上并不能被视为一个常数，因此该方法的除雾性能不稳定。LIU等^[13]提出基于多尺度相关小波方法，利用该方法在低频部分去雾，高频部分去噪并增强纹理细节，其能显著提高雾霾场景的感知可视性。GALDRAN等^[14]通过伽马校正从单一图像中提取多曝光图像，然后利用多尺度拉普拉斯混合方案将曝光后的图像融合为清晰图像。ZHENG等^[15]通过伽马校正和饱和度空间线性调整提取曝光不足的图像序列，利用基于自适应结构分解的MEF算法进行融合，从而提高局部细节的视觉质量。

现有方法大多基于深度信息估计进行去雾，去雾性能在很大程度上取决于精确的深度信息，为了实现高质量的去雾，需要在图像去雾过程中准确估计深度信息。FANG等^[16]提出一种全变分（Total Variation，TV）正则化变分模型，以同时估计透射图和去雾后图像，TV规则项能够较好地保持图像的边缘，利用TV项细化传输能够提高模糊图像的能见度。但是，该模型由于忽略了图像不同通道间的耦合，导致彩色图像的边缘模糊。WANG等^[17]针对边缘模糊的现象，采用3种不同耦合形式的变分模型（LTV、MTV、CTV）与暗原色先验相结合，有效改善了图像边缘模糊的现象，但是，基于一阶TV规则项的方法往往在边缘处产生阶梯伪影。高珠珠等^[18]利用二阶变分模型代替一阶TV项来约束透射图和去雾图像的估计，有效抑制了边缘阶梯伪影。

本文将基于TV规则项和BH（Bounded Hessian）规则项（以下简称TVBH）的变分模型与暗原色先验相结合，提出一种变分去雾模型H-TVBH。通过暗原色先验和四叉树分解估计初始透射率图和大气光值，将其代入H-TVBH模型中，结合分裂Bregman算法和快速傅里叶变换求得优化后的透射率和去雾图像。

1 相关工作 1.1 暗原色先验原理

在计算机视觉和图像处理中，大气散射模型^[19]被广泛应用于雾天图像形成的描述过程，具体如下：

$ I\left(x\right)=J\left(x\right)t\left(x\right)+A(1-t(x\left)\right) $

(1)

其中：$ I\left(x\right) $为输入有雾图像；$ J\left(x\right) $为恢复的清晰图像；$ A $为全局大气光值；$ t\left(x\right)(0\le t(x)\le 1) $为场景光的透射率。基于大气散射模型的去雾是从观测图像$ I\left(x\right) $中恢复原始无雾图像$ J\left(x\right) $，已知量仅为$ I\left(x\right) $，要求得$ J\left(x\right) $，需要透射率$ t\left(x\right) $和大气光值$ A $。

根据暗原色先验理论，在大多数户外无雾图像的非天空区域都存在亮度较低（接近于0）的通道。对于图像$ J $，暗通道表示为：

$ {J}^{\mathrm{d}\mathrm{a}\mathrm{r}\mathrm{k}}\left(x\right)=\underset{c\in \{\mathrm{R}, \mathrm{G}, \mathrm{B}\}}{\mathrm{m}\mathrm{i}\mathrm{n}}\left(\underset{y\in \mathrm{\Omega }\left(x\right)}{\mathrm{m}\mathrm{i}\mathrm{n}}\left({J}^{c}\right(y\left)\right)\right) $

(2)

其中：$ {J}^{\mathrm{d}\mathrm{a}\mathrm{r}\mathrm{k}} $代表室外无雾图像的暗通道；c代表R、G、B三通道；$ {J}^{c}\left(y\right) $代表恢复清晰图像的c通道图像；Ω（x）代表以x像素为中心的局部区域。根据观察到的统计数据，$ {J}^{\mathrm{d}\mathrm{a}\mathrm{r}\mathrm{k}} $的强度值始终很低并且接近于0。

将式（1）两边同时除以A，计算3个颜色通道的最小值，得到如下表达式：

$ \underset{y\in \mathrm{\Omega }\left(x\right)}{\mathrm{m}\mathrm{i}\mathrm{n}}\left(\underset{c}{\mathrm{m}\mathrm{i}\mathrm{n}}\frac{{I}^{c}\left(y\right)}{{A}^{c}}\right)=\tilde{t}\left(x\right)\underset{y\in \mathrm{\Omega }\left(x\right)}{\mathrm{m}\mathrm{i}\mathrm{n}}\left(\underset{c}{\mathrm{m}\mathrm{i}\mathrm{n}}\frac{{J}^{c}\left(y\right)}{{A}^{c}}\right)+1-\tilde{t}\left(x\right) $

(3)

结合暗原色先验可以估计出初始透射率$ \tilde{t}\left(x\right) $：

$ \tilde{t}\left(x\right)=1-\underset{y\in \mathrm{\Omega }\left(x\right)}{\mathrm{m}\mathrm{i}\mathrm{n}}\left(\underset{c}{\mathrm{m}\mathrm{i}\mathrm{n}}\frac{{I}^{c}\left(y\right)}{{A}^{c}}\right) $

(4)

为了使处理后的图像具有层次感和真实感，且保持图像的主观视觉效果，在式（4）中加入一个恒定的参数$ \omega (0 < \omega \le 1) $，该值根据实际情况而定，通常为0.90~0.97。

取雾天暗通道图像中的前0.1%像素，找到像素在原始图像中的对应位置，将具有最高亮度的点的值作为大气光值。根据计算得到的透射率$ t\left(x\right) $和大气光值$ A $，得到无雾图像$ J\left(x\right) $：

$ J\left(x\right)=\frac{I\left(x\right)-A}{\mathrm{m}\mathrm{a}\mathrm{x}\left(t\right(x), {t}_{0})}+A $

(5)

为了避免透射率t（x）趋近于0时图像颜色过饱和，可设置一个下限$ {t}_{0} $，通常令$ {t}_{0}=0.1 $。

1.2 TVBH模型

一阶TV规则项可以有效去除图像中的噪声，并保持图像边缘，但易产生边缘处阶梯效应和图像对比度低的问题；而二阶BH规则项在保持图像边缘的同时具有较好的平滑效果，抑制了阶梯效应。TVBH模型将上述两者相结合，可以较好地恢复图像的边缘细节，提高图像对比度，其模型形式如下：

$ E\left(u\right)=\alpha {\int }_{\mathrm{\Omega }}\left|\nabla u\right|\mathrm{d}x+\beta {\int }_{\mathrm{\Omega }}\left|{\nabla }^{2}u\right|\mathrm{d}x+\frac{1}{2}{\int }_{\mathrm{\Omega }}{\left(u-f\right)}^{2}\mathrm{d}x $

(6)

对于彩色图像，由于各层边缘扩散强度不同，易造成图像边缘模糊。为了使边缘信息得到较好的保持，可采用不同层图像之间耦合的形式得到能量泛函，如下：

$ \begin{array}{*{20}{c}} {E\left( u \right) = \alpha {\smallint _\Omega }\sqrt {\sum\limits_{i = 1}^3 {{{\left| {\nabla {u_i}} \right|}^2}} } {\rm{d}}x + \beta {\smallint _\Omega }\sqrt {\sum\limits_{i = 1}^3 {{{\left| {{\nabla ^2}{u_i}} \right|}^2}} } {\rm{d}}x + }\\ {\frac{1}{2}\sum\limits_{i = 1}^3 {{\smallint _\Omega }} {{\left( {{u_i} - {f_i}} \right)}^2}{\rm{d}}x} \end{array} $

(7)

其中：式（7）右侧前2项为TVBH模型的规则项，用以保持图像边缘细节；α、β为惩罚参数，决定处理后图像的光滑程度；最后一项为数据项，使恢复后的图像更接近于原始图像。

$ \nabla {u}_{i} $、$ {\nabla }^{2}{u}_{i} $定义分别为：

$ \nabla {u}_{i}=\left(\begin{array}{cc}{\partial }_{x}{u}_{i}& {\partial }_{y}{u}_{i}\end{array}\right) $

$ {\nabla }^{2}{u}_{i}=\left(\begin{array}{cc}{\partial }_{x}{\partial }_{x}{u}_{i}& {\partial }_{y}{\partial }_{x}{u}_{i}\\ {\partial }_{x}{\partial }_{y}{u}_{i}& {\partial }_{x}{\partial }_{y}{u}_{i}\end{array}\right) $

2 本文算法

为了抑制图像中的噪声并充分保留去雾后图像的边缘特征，本文将TVBH模型与暗原色先验相结合，从而实现图像去雾。本文新的能量泛函形式为：

$ \begin{array}{*{20}{c}} {E(t,u) = \lambda {\smallint _\Omega }\left| {\nabla t} \right|{\rm{d}}x + \alpha {\smallint _\Omega }\sqrt {\sum\limits_{i = 1}^3 {{{\left| {\nabla {u_i}} \right|}^2}} } {\rm{d}}x + }\\ {\frac{1}{2}{\smallint _\Omega }{{\left( {t - \tilde t} \right)}^2}{\rm{d}}x + \beta {\smallint _\Omega }\sqrt {\sum\limits_{i = 1}^3 {{{\left| {{\nabla ^2}{u_i}} \right|}^2}} } {\rm{d}}x + }\\ {\frac{1}{2}\sum\limits_{i = 1}^3 {{\smallint _\Omega }} {{\left( {{f_i} - {u_i} \cdot t - (1 - t)A} \right)}^2}{\rm{d}}x} \end{array} $

(8)

其中：t为待求解的透射率；$ \tilde{t} $为初始透射率图；$ A $为大气光值；$ {\int }_{\mathrm{\Omega }}\left|\nabla t\right|\mathrm{d}x $是透射率的规则项，$ \frac{1}{2}{\int }_{\mathrm{\Omega }}{\left(t-\tilde{t}\right)}^{2}\mathrm{d}x $是透射率的数据项，通过最小化这2项来建模优化透射率；λ为透射率规则项的惩罚参数，用以保证透射率图像的光滑程度；$ {\smallint _\Omega }\sqrt {\sum\limits_{i = 1}^3 {{{\left| {\nabla {u_i}} \right|}^2}} } {\rm{d}}x $和$ {\smallint _\Omega }\sqrt {\sum\limits_{i = 1}^3 {{{\left| {{\nabla ^2}{u_i}} \right|}^2}} } {\rm{d}}x $是去雾图像的规则项，用以保持恢复后图像的边缘；$ \sum\limits_{i = 1}^3 {{\smallint _\Omega }} {\left( {{f_i} - {u_i} \cdot t - (1 - t)A} \right)^2}{\rm{d}}x $为去雾图像的数据项，使去雾后图像更接近于原图像；α、β为去雾图像规则项的惩罚参数，用以保证去雾后图像的光滑程度。

2.1 大气光值估计

通常情况下将图像中最亮颜色的值估计为大气光值，然而，在对图像进行大气光值估计时，会由于图像中存在过亮物体或天空区域而导致大气光值估计错误。因此，本文采用四叉树分解法^[20]，将原图像分成4个均匀的子区块，对子区块再进行分割，分割过程是迭代而不断重复的，直到满足一定的阈值标准（阈值标准为长×宽=200像素的区域）。对每个区块的大气光值使用暗原色先验进行计算，即将最小滤波器应用到每个块，然后将最大值作为图像的全局大气光值。

2.2 模型求解

为了提高模型的计算效率，本文采用基于快速傅里叶变换^[21]的分裂Bregman算法，引入辅助变量$ \boldsymbol{x} $，$ \boldsymbol{w}=({\boldsymbol{w}}_{1}, {\boldsymbol{w}}_{2}, {\boldsymbol{w}}_{3}{)}^{\mathrm{T}} $，$ \boldsymbol{v}=({\boldsymbol{v}}_{1}, {\boldsymbol{v}}_{2}, {\boldsymbol{v}}_{3}{)}^{\mathrm{T}} $，令$ \boldsymbol{x}=\nabla t $，$ \boldsymbol{w}=\nabla u $，$ \boldsymbol{v}={\nabla }^{2}u $，则上述能量泛函转化为：

$ \begin{array}{*{20}{l}} {E(t,{\boldsymbol{x}},u,{\boldsymbol{w}},{\boldsymbol{v}}) = }\\ {\lambda {\smallint _\Omega }\left| {\boldsymbol{x}} \right|{\rm{d}}x + \frac{1}{2}{\smallint _\Omega }{{\left( {t - \tilde t} \right)}^2}{\rm{d}}x + }\\ {\frac{1}{2}\sum\limits_{i = 1}^3 {{\smallint _\Omega }} {{\left( {{f_i} - {u_i} \cdot t - (1 - t)A} \right)}^2}{\rm{d}}x + }\\ {\frac{{{\mu _1}}}{2}{\smallint _\Omega }{{({\boldsymbol{x}} - \nabla t - {\boldsymbol{q}})}^2}{\rm{d}}x + \alpha {\smallint _\Omega }\sqrt {\sum\limits_{i = 1}^3 {{{\left| {{{\boldsymbol{w}}_{\boldsymbol{i}}}} \right|}^2}} } {\rm{d}}x + }\\ {\beta {\smallint _\Omega }\sqrt {\sum\limits_{i = 1}^3 {{{\left| {{{\boldsymbol{v}}_{\boldsymbol{i}}}} \right|}^2}} } {\rm{d}}x + \frac{{{\mu _2}}}{2}\sum\limits_{i = 1}^3 {{\smallint _\Omega }} {{({{\boldsymbol{w}}_{\boldsymbol{i}}} - \nabla {u_i} - {{\boldsymbol{b}}_{\boldsymbol{i}}})}^2}{\rm{d}}x + }\\ {\frac{{{\mu _3}}}{2}\sum\limits_{i = 1}^3 {{\smallint _\Omega }} {{({{\boldsymbol{v}}_i} - {\nabla ^2}{u_i} - {{\boldsymbol{d}}_i})}^2}{\rm{d}}x} \end{array} $

(9)

其中：μ₁、μ₂、μ₃为惩罚参数；$ \boldsymbol{q}\boldsymbol{、}{\boldsymbol{b}}_{i}\mathrm{、}{\boldsymbol{d}}_{i} $为Bregman参数。

迭代方式为：

$ \left\{\begin{array}{l}{\boldsymbol{q}}_{}^{k+1}={\boldsymbol{q}}_{}^{k}+\nabla t-\boldsymbol{x}\\ {\boldsymbol{b}}_{i}^{k+1}={\boldsymbol{b}}_{i}^{k}+\nabla {u}_{i}-{\boldsymbol{w}}_{i}\\ {\boldsymbol{d}}_{i}^{k+1}={\boldsymbol{d}}_{i}^{k}+{\nabla }^{2}{u}_{i}-{\boldsymbol{v}}_{i}\end{array}\right. $

(10)

对上述能量泛函交替优化求极值，固定其中4个变量求另一个变量，交替迭代求t、$ \boldsymbol{x} $、$ u $、$ \boldsymbol{v} $、$ \boldsymbol{w} $的最小值，则子问题的极值问题转化为：固定$ {\boldsymbol{x}}^{k}, {u}^{k}, {\boldsymbol{w}}^{k}, {\boldsymbol{v}}^{k} $，求t；固定$ {t}^{k+1}, {u}^{k}, {\boldsymbol{w}}^{k}, {\boldsymbol{v}}^{k} $，求$ \boldsymbol{x} $；固定$ {t}^{k+1}, {\boldsymbol{x}}^{k+1}, {\boldsymbol{w}}^{k}, {\boldsymbol{v}}^{k} $，求u；固定$ {t}^{k+1}, {\boldsymbol{x}}^{k+1}, {u}^{k+1}, {\boldsymbol{v}}^{k} $，求$ \boldsymbol{w} $；固定$ {t}^{k+1}, {\boldsymbol{x}}^{k+1}, {\boldsymbol{u}}^{k+1}, $$ {\boldsymbol{w}}^{\boldsymbol{k}+1} $，求$ \boldsymbol{v} $。由此得到关于t、$ \boldsymbol{x} $、u、$ \boldsymbol{v} $、$ \boldsymbol{w} $的表达式如下：

$ \begin{array}{l}t+t(A-{u}^{k}{)}^{2}-{\mu }_{1}\nabla \cdot \left(\nabla t\right)=\\ \tilde{t}-{\mu }_{1}\nabla \cdot (\boldsymbol{x}-\boldsymbol{q})-(f-A)(A-{u}^{k})\end{array} $

(11)

$ \boldsymbol{x}=\mathrm{m}\mathrm{a}\mathrm{x}\left\{\left|\nabla {t}^{k+1}+\boldsymbol{q}\right|-\frac{\lambda }{{\mu }_{1}}, 0\right\}\cdot \frac{\nabla {t}^{k+1}+\boldsymbol{q}}{\left|\nabla {t}^{k+1}+\boldsymbol{q}\right|} $

(12)

$ \begin{array}{l}{u}_{i}\cdot {t}^{2}-{\mu }_{2}\nabla \cdot \left(\nabla {u}_{i}\right)+{\mu }_{3}{\nabla }^{2}\cdot \left({\nabla }^{2}{u}_{i}\right)=\\ {f}_{i}\cdot t-A(1-t)t-{\mu }_{2}\nabla \cdot \left({\boldsymbol{w}}_{\boldsymbol{i}}-{\boldsymbol{b}}_{\boldsymbol{i}}\right)+{\mu }_{3}{\nabla }^{2}\cdot \left({\boldsymbol{v}}_{\boldsymbol{i}}-{\boldsymbol{d}}_{\boldsymbol{i}}\right)\end{array} $

(13)

$ {\boldsymbol{w}}_{i}^{k+1}=\mathrm{m}\mathrm{a}\mathrm{x}\left\{\left|\nabla {u}_{i}^{k+1}+{\boldsymbol{b}}_{i}\right|-\frac{\alpha }{{\mu }_{2}}, 0\right\}\cdot \frac{\nabla {u}_{i}^{k+1}+{\boldsymbol{b}}_{i}}{\left|\nabla {u}_{i}^{k+1}+{\boldsymbol{b}}_{i}\right|} $

(14)

$ {\boldsymbol{v}}_{i}^{k+1}=\mathrm{m}\mathrm{a}\mathrm{x}\left\{\left|{\nabla }^{2}{u}_{i}^{k+1}+{\boldsymbol{d}}_{\boldsymbol{i}}\right|-\frac{\beta }{{\mu }_{3}}, 0\right\}\cdot \frac{{\nabla }^{2}{u}_{i}^{k+1}+{\boldsymbol{d}}_{i}}{\left|{\nabla }^{2}{u}_{i}^{k+1}+{\boldsymbol{d}}_{i}\right|} $

(15)

其中：式（11）、式（13）采用离散傅里叶变换进一步优化，从而提高算法的计算效率；式（12）、式（14）、式（15）为广义软阈值公式，$ \left| {\nabla u_i^{k + 1} + {{\boldsymbol{b}}_i}} \right| = \sqrt {\sum\limits_{i = 1}^3 {{{\left| {\nabla u_i^{k + 1} + {{\boldsymbol{b}}_i}} \right|}^2}} } $，$ \left| {{\nabla ^2}u_i^{k + 1} + {{\boldsymbol{d}}_i}} \right| = \sqrt {\sum\limits_{i = 1}^3 {{{\left| {{\nabla ^2}u_i^{k + 1} + {{\boldsymbol{d}}_i}} \right|}^2}} } $，基于耦合的方式使边缘保持效果更好。

对式（11）、式（13）进行傅里叶变换，得到t和u的解析解分别如下：

$ {t}_{i, j}=\mathfrak{R}\left(\mathrm{F}\mathrm{F}{\mathrm{T}}^{-1}\left(\frac{\mathrm{F}\mathrm{F}\mathrm{T}\left({G}_{i, j}\right)}{{\xi }_{1}}\right)\right) $

(16)

$ {u}_{i, j}=\mathfrak{R}\left(\mathrm{F}\mathrm{F}{\mathrm{T}}^{-1}\left(\frac{\mathrm{F}\mathrm{F}\mathrm{T}\left({H}_{i, j}\right)}{{\xi }_{2}}\right)\right) $

(17)

其中：$ {\xi }_{1} $、$ {\xi }_{2} $是等式左边通过离散傅里叶变换后基于频率域的等价变换，分别为$ {\xi }_{1}=1+(A-{u}^{k}{)}^{2}-2{\mu }_{1}\left(\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{2\mathrm{\pi }s}{N}\right)+\right. $$ \left.\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{2\mathrm{\pi }r}{M}\right)-2\right) $，$ {\xi }_{2}={t}^{2}-2{\mu }_{2}\left(\mathrm{c}\mathrm{o}\mathrm{s}\frac{2\mathrm{\pi }s}{N}+\right. $ $ \left.\mathrm{c}\mathrm{o}\mathrm{s}\frac{2\mathrm{\pi }r}{M}-2\right)+4{\mu }_{3}{\left(\mathrm{c}\mathrm{o}\mathrm{s}\frac{2\mathrm{\pi }s}{N}+\mathrm{c}\mathrm{o}\mathrm{s}\frac{2\mathrm{\pi }r}{M}-2\right)}^{2} $；$ {G}_{i, j} $和$ {H}_{i, j} $为等式右边的离散形式，分别为$ {G}_{i, j}=\tilde{t}-{\mu }_{1}\left({\partial }_{x}^{-}\right({\boldsymbol{x}}_{1i, j}-{\boldsymbol{q}}_{1i, j})+{\partial }_{y}^{-}({\boldsymbol{x}}_{2i, j}-{\boldsymbol{q}}_{2i, j}\left)\right)-(f-A)(A-{u}^{k}) $，$ {H_{i,j}} = t - A(1 - t)t - {\mu _2}(\partial _x^ - ({{\boldsymbol{w}}_{1i,j}} - {{\boldsymbol{b}}_{1i,j}}) + \partial _y^ - ({{\boldsymbol{w}}_{2i,j}} - {{\boldsymbol{b}}_{2i,j}})) + {\mu _3}(\partial _x^ + \partial _x^ - ({{\boldsymbol{v}}_{1i,j}} - {{\boldsymbol{d}}_{1i,j}}) $ $ + \partial _y^ - \partial _x^ - ({{\boldsymbol{v}}_{2i,j}} - {{\boldsymbol{d}}_{2i,j}}) + \partial _x^ - \partial _y^ - ({{\boldsymbol{v}}_{3i,j}} - {{\boldsymbol{d}}_{3i,j}}) + \partial _y^ + \partial _y^ - ({{\boldsymbol{v}}_{4i,j}} - {{\boldsymbol{d}}_{4i,j}}))$；i、r、j、s为索引，$ i, r=\mathrm{0, 1}, \cdot \cdot \cdot , M-1;j, s=\mathrm{0, 1}, \cdots , N-1 $。

本文算法流程如图 1所示。

本文算法输入为原始图像f，输出为去雾图像u，求解步骤如下：1）根据暗原色先验估计图像原始透射率$ \tilde{t} $；2）利用四叉树分解法估计大气光值A；3）初始化，将$ \boldsymbol{x}\boldsymbol{、}{\boldsymbol{w}}_{\boldsymbol{i}}\boldsymbol{、}{\boldsymbol{v}}_{\boldsymbol{i}}\boldsymbol{、}\boldsymbol{q}\boldsymbol{、}{\boldsymbol{b}}_{\boldsymbol{i}}\boldsymbol{、}{\boldsymbol{d}}_{\boldsymbol{i}} $值设为0，惩罚参数$ \lambda \boldsymbol{、}\alpha \boldsymbol{、}\beta \boldsymbol{、}{\mu }_{1}\boldsymbol{、} $ $ {\mu }_{2}\boldsymbol{、}{\mu }_{3} > 0 $；4）根据$ \tilde{t} $、A和原始图像f，迭代求解各变量；5）收敛直至得到t和u。

3 实验结果与分析

本文实验在Windows10系统、Intel Core i5-6500处理器、8 GB RAM、使用Matlab R2014a的PC机上完成。为了验证所提方法的有效性，从LIVE Image Defogging数据库^[22]中选取几幅图像模拟有雾图像，将本文方法与其他几种经典去雾算法进行定性和定量比较。实验主要参数设置为：$ \omega =0.95 $，$ \lambda =0.01 $，$ \alpha =0.25 $，$ \beta =0.5 $，$ {\mu }_{1}=0.05 $，$ {\mu }_{2}=0.000\mathrm{ }1 $，$ {\mu }_{3}=0.000\mathrm{ }005 $。其中，$ \lambda $和$ {\mu }_{1} $保证透射率图的光滑程度，$ \alpha $和$ \beta $保证去雾后图像的光滑程度，$ {\mu }_{2} $、$ {\mu }_{3} $值设置较小，若设置过大会导致图像模糊并丢失细节特征。

3.1 主观评价

图 2所示为各方法对航空图、小麦图、森林图、天空图以及局部放大图的去雾效果对比。从图 2可以看出：对于航空图及放大区域，各方法都可以在不同程度上去雾，但文献[7]方法、文献[13]方法、文献[14]方法在远景处仍存在残留的雾，文献[14]方法过于平滑，在一定程度上模糊了纹理，文献[5]方法近景部分对比度相对较低，本文方法无论是在近景还是远景处均更加清晰地恢复了场景细节；对于小麦图，文献[5]方法、文献[13]方法去雾后图像对比度和饱和度较低，文献[7]方法、文献[14]方法提高了图像亮度，但远景部分（方框处）有明显的雾残留，本文方法增强了远景可见性，具有更好的去雾效果；对于森林图及放大区域，各方法均有一定的雾残留，文献[5]方法去雾后图像颜色过饱和，出现颜色失真的现象，文献[13]方法去雾力度明显不足，有大量雾残留，本文方法相对其他4种方法去雾效果较好；在天空图中，文献[5]方法天空区域较为灰暗，文献[7]方法、文献[13]方法提高了图像对比度和亮度，但部分纹理细节丢失，文献[14]方法去雾后图像出现光晕伪影，本文方法保留了更多的纹理细节，且增强了远景可见性，视觉效果更为自然。图 3所示为更多有雾图像的去雾效果，从图 3可以看出：文献[5]方法去雾后图像过饱和，颜色失真，且近景区域（方框处）去雾不彻底；文献[7]方法、文献[14]方法图像对比度有所提高，但远处仍有一些残留雾，文献[14]方法去雾后图像细节有些模糊；文献[13]方法处理后的图像整体亮度较为灰暗，边缘细节丢失严重，视觉效果不自然；本文方法去雾后的图像有效保留了边缘细节，提高了图像对比度和饱和度，但部分图像颜色偏暗。

图 4所示为各方法对玩偶有雾图像和条纹有雾图像放大区域的去雾效果对比，将本文方法与文献[16-17]中的LTV、MTV、CTV这3种方法以及文献[18]中的H-TGV方法，针对图像边缘以及噪声去除性能进行比较。由图 4可以看出：对于玩偶图，文献[16]方法提高了图像对比度，但去雾后噪声放大，LTV^[17]方法去雾后图像纹理细节模糊不清，且图像颜色过饱和，MTV^[17]和CTV^[17]2种方法去雾效果相对较好，但图像略平滑，丢失部分边缘细节，H-TGV^[18]和本文方法去雾效果接近，均有效改善了图像的饱和度和对比度，保留了图像的纹理细节；对于含噪声的条纹图，文献[16]方法和LTV^[17]方法处理后的图像仍有噪声残留，MTV^[17]和CTV^[17]方法去噪效果相对较好，但中间边缘处模糊不清，H-TGV^[18]和本文方法在去雾的同时有效抑制了噪声，且边缘处细节保持较好，去雾效果更加突出。

3.2 客观评价

为了更好地验证本文方法的有效性和鲁棒性，采用新增可见边之比（e）、可见边规范化梯度之比（$ \stackrel{-}{r} $）、雾密度感知（D）、结构相似性（SSIM）、峰值信噪比（PSNR）作为评价指标。对图 2、图 3采用新增可见边e、可见边规范化之比$ \stackrel{-}{r} $和雾密度感知D进行评价，其中，e、$ \stackrel{-}{r} $值越大，D值越小，表明图像去雾效果越好。表 1、表 2所示为图 2、图 3各方法的去雾结果对比。从表 1、表 2可以看出，文献[5]方法、文献[13]方法中的$ \stackrel{-}{r} $结果偏大，这是因为去雾后图像过增强，颜色过饱和，而本文方法具有较大的e和$ \stackrel{-}{r} $以及较小的雾密度感知D，表明本文方法具有较好的去雾性能。

对于模拟有雾图像，由于已有相应的清晰图像作参考，因此可以采用SSIM和PSNR，通过不同的策略来进一步评价结果的准确性。通常情况下，SSIM值越大，表示去雾结果与真实图像的结构相似性越好，去雾图像越接近于真实图像；PSNR值越大，表明图像去雾后的效果更加理想。此外，本文还对不同变分去雾模型的运行时间进行了比较。表 3所示为各方法对图 4的去雾效果对比。从表 3可以看出：本文方法可以获得较高的PSNR值和SSIM值，去雾效果较好，去雾图像更接近于真实图像；虽然H-TGV^[18]方法在PSNR和SSIM值上与本文方法接近，但运行时间相对较长；与其他方法相比，本文方法具有较短的运行时间。

4 结束语

本文提出一种结合暗原色先验与TVBH变分模型的去雾方法，根据暗原色先验估计透射率图，通过四叉树分解得到大气光值，利用分裂Bregman算法和快速傅里叶变换求得优化后的透射率和去雾图像。实验结果表明，该方法能够较好地保留图像的边缘细节，改善图像的可视化效果。但是，本文方法去雾后图像颜色偏暗，解决该问题并提高去雾图像的视觉质量将是下一步的研究方向。