Single Image Haze Removal Using Dark Channel Prior

工作简述

何神的第一篇经典。具体工作如下图所示。开创暗通道先验来做图像去雾的工作。直接进入方法解析部分。

方法论

计算机大气透射模型

在计算机视觉和计算机图形学中，广泛用于描述雾霾图像形成的模型如下：

其中是观测到的强度，是场景辐射亮度，是全局大气光，是描述未散射并到达相机的光部分的介质透射率。

去雾的目标是从中恢复、和。方程（1）右侧的第一项称为直接衰减 [16]，第二项称为空气光 [6, 16]。直接衰减描述了场景辐射及其在介质中的衰减，而空气光则由先前散射的光产生并导致场景颜色的变化。当大气均匀时，透射率t可表示为：

其中是大气的散射系数。它表明场景辐射亮度随着场景深度呈指数衰减。

简而言之，目标是从中恢复。而根据这个方程而言，我们需要估算和。

暗通道先验

这个网页没做min的下标渲染导致只能截图

其中是的颜色通道，是以为中心的局部补丁。我们的观察表明，如果是无雾的室外图像，则除天空区域外，的强度较低且趋于零。我们将称为的暗通道，并将上述统计观察或知识称为暗通道先验.

暗通道中的低强度主要是由于三个因素造成的：

a) 阴影。例如，城市景观图像中汽车、建筑物和窗户内部的阴影，或者风景图像中树叶、树木和岩石的阴影；

b) 彩色物体或表面。例如，任何颜色通道中缺乏颜色的任何物体（例如，绿草/树木/植物、红色或黄色花/叶以及蓝色水面）将导致暗通道中的值较低；

c) 深色物体或表面。例如，深色的树干和石头。由于自然的户外图像通常充满阴影且色彩丰富，因此这些图像的暗通道非常暗！

即，对于一般的正常图像而言，rgb空间总有一个通道有较低的值（且趋近为0）（除了白色的天空）。根据这个来估算暗通道先验，进一步估算大气的透过率。

根据暗通道先验还原图像

这里，首先假设给定大气光。后面再额外介绍一种自动估计大气光的方法。我们进一步假设局部贴片中的传输是恒定的。我们将补丁的传输表示为。对前置的雾霾成像方程进行局部补丁中的min运算，我们有：

公式推导组合技

这样得到了传输补丁的估计。然后对于前面的天空区域，虽然不适用暗通道先验，但是天空成分和大气散射光一致，即天空的颜色通常与雾霾图像中的大气光非常相似。由于天空无限大并且透射率往往为零，因此上面的等式（11）可以很好地处理天空区域和非天空区域。不需要事先分离天空区域。

实际上，即使在晴天，大气中也并非完全没有任何颗粒。所以，当我们看远处的物体时，雾霾仍然存在。此外，雾霾的存在是人类感知深度的基本线索[3, 13]。这种现象称为空中透视（怎么和画画时候要考虑的一样）。如果彻底去除雾气，图像可能会显得不自然，并且可能会失去深度感。因此，我们可以通过在方程（11）中引入一个常数参数 ω (0<ω≤1) 来选择为远处的物体保留非常少量的雾度：

它相当不错，但包含一些块效应，因为补丁中的透射率并不总是恒定的。(也就是透射率的图比较粗糙) 在下一小节中，我们使用soft matting方法细化该贴图。

而对这个估算的解释也很简单。即在一个区域（卷积核）内的像素的颜色通道中最低的颜色空间除以大气散射光在被1相减即得到估算的大气透过率t。

soft matting

这里的b图像可以看出估计的传输图是相当粗糙的，想得到更加精细的去雾图像，则需要进一步处理。

该方法在前置的方法上进一步修正得到更好的估算t，其缺点是速度特别慢，不适用在实时场合，2011年，又提出可以使用导向滤波的方式来获得更细腻的结果，这个方法的运算主要集中在方框滤波（均值滤波），而这种操作在opencv或者其他的图像库中都有快速算法。可以考虑使用

首先，把精细化任务改为最小化如下误差方程：

其中是Levin[7]提出的Matting Laplacian矩阵，是正则化参数。第一项是平滑项，第二项是数据项。

后续描述如上图所示，不赘述思想，是Levin的软抠图方法。

图像恢复

得到了精细的大气透过率图片之后，最后根据公式，得到下面的图像恢复方程。

其中，的典型值是0.1

大气光的估计

关于A的估计，过去的方法主要是取图像当中最亮的点。当然这样会带来一个问题：被白色的建筑物和一些白色汽车的反光给标记上。

而何神继续使用暗通道来改进大气光估计：

首先选取暗通道中最亮的 0.1% 像素。这些像素是最不透明的。在这些像素中，选择输入图像I中强度最高的像素作为大气光。

如图中的红色框框就是选取出来暗通道图像当中最亮的部分。基本上可以定位到大气当中。

实验讨论

对于何神的分析，这边还需要进一步的实验讨论。复现其代码看看效果：

参考博客：暗通道去雾算法原理及实现

首先参考这个博客内的第二个代码，对每一步有相当详细的函数封装：(这边我再加入一些注释)——

import cv2
import math
import numpy as np

# 暗通道先验（输入原图像和卷积核的大小）
def DarkChannel(im,sz):
    b,g,r = cv2.split(im)          # 颜色通道拆分
    dc = cv2.min(cv2.min(r,g),b)   # 寻找通道中的最小值
    # 设置一个矩形卷积核
    kernel = cv2.getStructuringElement(cv2.MORPH_RECT,(sz,sz))
    dark = cv2.erode(dc,kernel)    # 腐蚀，（类似最小值滤波的过程（即一个补丁内取最小的rgb值））
    return dark                    # 返回暗通道图片（单通道的图）

# 估计大气光（输入原图像）
def AtmLight(im,dark):
    [h,w] = im.shape[:2]
    imsz = h*w                      # 求全局像素个数（图像大小）
    numpx = int(max(math.floor(imsz/1000),1))    # 前千分之一的暗通道像素数量（不低于1）
    darkvec = dark.reshape(imsz)    # 暗通道图片拉直为向量
    imvec = im.reshape(imsz,3)      # 原图像转为3颜色通道的，像素拉直
  
    indices = darkvec.argsort()     # 返回暗通道从小到大的索引值
    indices = indices[imsz-numpx::] # 选择暗通道图像值最大（最亮，透射率最低）的部分
  
    atmsum = np.zeros([1,3])        
    for ind in range(1,numpx):
       atmsum = atmsum + imvec[indices[ind]] # 前千分之一的像素值求和
  
    A = atmsum / numpx # 取平均值，得到大气光
    return A

# 得到大气透过率估计值（纯公式，见上面）
def TransmissionEstimate(im,A,sz):
    omega = 0.95
    im3 = np.empty(im.shape,im.dtype)
  
    for ind in range(0,3):
        im3[:,:,ind] = im[:,:,ind]/A[0,ind]
  
    transmission = 1 - omega*DarkChannel(im3,sz)
    return transmission

# 使用导向滤波的方式来获得更细腻的结果，这个方法的运算主要集中在方框滤波（均值滤波）
# 传入原灰度图，大气传输率图像（待滤波），窗口滤波核大小，极小值（防止除以0）
def Guidedfilter(im,p,r,eps):
    # 计算im,p的均值
    mean_I = cv2.boxFilter(im,cv2.CV_64F,(r,r))
    mean_p = cv2.boxFilter(p, cv2.CV_64F,(r,r))
    # 计算乘积的均值，进一步得到协方差
    mean_Ip = cv2.boxFilter(im*p,cv2.CV_64F,(r,r))
    cov_Ip = mean_Ip - mean_I*mean_p
    # 计算im的平方均值，进一步得到方差
    mean_II = cv2.boxFilter(im*im,cv2.CV_64F,(r,r))
    var_I   = mean_II - mean_I*mean_I
  
    # 使用协方差和方差计算线性系数a和b
    a = cov_Ip/(var_I + eps)
    b = mean_p - a*mean_I
    # 对a和b均值滤波
    mean_a = cv2.boxFilter(a,cv2.CV_64F,(r,r))
    mean_b = cv2.boxFilter(b,cv2.CV_64F,(r,r))
  
    return mean_a*im + mean_b

def TransmissionRefine(im,et):
    # 转化为灰度图
    gray = cv2.cvtColor(im,cv2.COLOR_BGR2GRAY)
    gray = np.float64(gray)/255
    r = 60
    eps = 0.0001
    # 得到精修的大气传输率图
    return Guidedfilter(gray,et,r,eps)

# tx 保留一部分大气模糊（体现图像的深度）
def Recover(im,t,A,tx = 0.1):
    res = np.empty(im.shape,im.dtype)
    t = cv2.max(t,tx)
  
    for ind in range(0,3):
        res[:,:,ind] = (im[:,:,ind]-A[0,ind])/t + A[0,ind]
  
    return res

if __name__ == '__main__':
  
    fn = './image/15.png'
  
    src = cv2.imread(fn)
    # 先化为[0,1]的空间（转化为小数方便运算）
    I = src.astype('float64')/255
  
    dark = DarkChannel(I,7)
    A = AtmLight(I,dark)
    te = TransmissionEstimate(I,A,15)
    # 下面的这个操作是对t的矫正（比较耗时）
    t = TransmissionRefine(src,te)
    J = Recover(I,t,A,0.1)
    cv2.imshow('J',J)
    cv2.waitKey()

关于实际的效果对比：

首先是原图和正常的处理（随便扒了一个网图）：

原图，扒了个网图

暗通道法去雾

然后按照其11年提出的方法做的t值的矫正删除之后：

纯暗通道先验

可以看见，暗通道先验主要起到的是一个对T的估计作用，而后续的soft matting起到关键的精修（当然这里选择的是其后续提出的更加高效的算法）

然后一整套流程处理的各个图片如下：

暗通道先验，大气透过率估计，精确估计

而耗时来说，11年用导向滤波改进的soft matting耗时不比暗通道估计耗时长（实际上是更短），在python下测试的情况是0.017s。是一种比较良好的方法（实时且效果良好）

大气散射模型

题外话，关于暗通道的方法，额外关注一下实际的大气光学散射公式，思考改进的空间。

1. 米氏散射模型

米氏散射模型用于描述光在大气中与中等大小颗粒（如雾霾中的水滴、灰尘等）发生的散射现象。该模型的散射强度依赖于激光的波长和颗粒的大小。

1.1 散射强度公式

米氏散射的散射强度通常由以下公式描述：

其中：

k = 是激光的波数，为激光波长（905 nm）。
是米氏散射的散射幅度。
和是散射角度。
是颗粒的散射相函数，依赖于颗粒的形状和光的传播方向。
是单位固体角内的散射强度。

1.2 计算步骤

计算光的传播方向和粒子散射的方向：确定入射激光与颗粒相互作用的角度。
选择适当的散射函数：根据粒子的大小、形状、折射率等选择合适的散射函数。
数值解算：由于米氏散射涉及到复杂的积分运算，通常使用数值方法（例如，有限差分法、有限元法）来求解散射强度。

1.3 相关计算工具

常用的计算工具包括MiePlot（专门用于米氏散射计算）或MATLAB中的散射模型库。

2. 吸收模型（Beer-Lambert定律）

吸收模型基于Beer-Lambert定律，用于描述激光光束在大气中传播时的衰减情况。

2.1 Beer-Lambert定律

Beer-Lambert定律的基本形式为：其中：

是距离源点为的位置处的光强。
是初始光强。
是吸收系数，描述了大气中材料对激光的吸收程度。它通常与大气的粒子浓度和激光波长有关。
是光束传播的距离。

2.2 吸收系数的计算

吸收系数可通过以下方式近似：其中：

是激光与颗粒或分子相互作用的有效吸收截面（与粒子的种类和大小有关）。
是大气中颗粒的浓度。

2.3 计算步骤

测定吸收系数：基于大气中颗粒的类型（如水蒸气、烟雾颗粒等）和浓度来估算吸收系数。
根据传播距离进行衰减计算：利用Beer-Lambert定律计算光强在不同传播距离下的变化。

2.4 相关计算工具

吸收系数的计算可以使用MODTRAN（Moderate Resolution Atmospheric Transmission）等大气传输模型工具来进行。
对于更简单的场景，也可以通过文献中的标准吸收系数数据来进行近似计算。

总结与解决方法

米氏散射模型：适用于雾霾中粒子大小与激光波长接近的情况，通常通过数值方法计算散射强度。可以借助现有的工具（如MiePlot）来进行计算。
吸收模型：利用Beer-Lambert定律描述光强随着传播距离的衰减。在复杂场景下，吸收系数可以通过实验或大气模型（如MODTRAN）获得。

去雾block：暗通道先验法和大气散射模型