1. 算法原理
基本原理来源于何凯明大神的CVPR09的论文Single Image Haze Removal Using Dark Channel Prior
暗通道
所谓暗通道是一个基本假设，这个假设认为，在绝大多数的非天空的局部区域中，某一些像素总会有至少一个颜色通道具有很低的值。这个其实很容易理解，实际生活中造成这个假设的原因有很多，比如汽车，建筑物或者城市中的阴影，或者说色彩鲜艳的物体或表面(比如绿色的树叶，各种鲜艳的花，或者蓝色绿色的睡眠)，颜色较暗的物体或者表面，这些景物的暗通道总是变现为比较暗的状态。
所以暗通道是什么呢？其实比较简单，作者认为暗通道是：
暗通道先验理论指出：
暗通道实际上是在rgb三个通道中取最小值组成灰度图，然后再进行一个最小值滤波得到的。我们来看一下有雾图像和无雾图像暗通道的区别：
可以发现，有雾的时候会呈现一定的灰色，而无雾的时候咋会呈现大量的黑色（像素为接近0），作者统计了5000多副图像的特征，基本都符合这样一条先验定理。
雾图形成模型
计算机视觉中，下面这个雾图形成模型是被广泛使用的：
其中I(x)是现有的图像（待去雾），J(x)是要恢复的原无雾图像，A是全球大气光成分，t(x)是透射率，现在的条件就是已知I(x)，来求J(x),显然不加任何限制的话是有无穷多个解的。
但是现实生活中，即使是晴天白云，空气中也会存在一些颗粒，看远方的物体还是能够感觉到雾的影响，另外，雾的存在可以让人们感觉到景深的存在，所以我们保留一部分的雾，上式修正为：其中w是[0-1]之间的一个值，一般取0.95差不多。
上面的推导都是假设全球大气光是已知的，实际中，我们可以借助暗通道图来从有雾图像中来获取该值：
1.从暗通道图中按照亮度大小取前0.1%的像素。
2.在这些位置中，在原始图像中寻找对应具有最高亮度点的值，作为A值。
到这里，我们就可以进行无雾图像的恢复了：
当投射图t很小时，会导致J的值偏大，会导致图片某些地方过爆，所以一般可以设置一个阈值来限制，我们设置一个阈值：一般设置较小，0.1即可。
利用这个理论的去雾效果就不错了，下面是我在网上找的例子：
但是这个去雾效果还是挺粗糙的，主要原因是由于透射率图过于粗糙了，何凯明在文章中提出了soft matting方法，然后其缺点是速度特别慢，不适用在实时场合，2011年，又提出可以使用导向滤波的方式来获得更细腻的结果，这个方法的运算主要集中在方框滤波（均值滤波），而这种操作在opencv或者其他的图像库中都有快速算法。可以考虑使用。
2.代码实现
# -*- coding: utf-8 -*-
"""
Created on Sat Jun 9 11:28:14 2018
@author: zhxing
"""
import cv2
import numpy as np
def zmMinFilterGray(src, r=7):
'''''最小值滤波，r是滤波器半径'''
return cv2.erode(src,np.ones((2*r-1,2*r-1)))
# =============================================================================
# if r <= 0:
# return src
# h, w = src.shape[:2]
# I = src
# res = np.minimum(I , I[[0]+range(h-1) , :])
# res = np.minimum(res, I[range(1,h)+[h-1], :])
# I = res
# res = np.minimum(I , I[:, [0]+range(w-1)])
# res = np.minimum(res, I[:, range(1,w)+[w-1]])
# =============================================================================
# return zmMinFilterGray(res, r-1)
def guidedfilter(I, p, r, eps):
'''''引导滤波，直接参考网上的matlab代码'''
height, width = I.shape
m_I = cv2.boxFilter(I, -1, (r,r))
m_p = cv2.boxFilter(p, -1, (r,r))
m_Ip = cv2.boxFilter(I*p, -1, (r,r))
cov_Ip = m_Ip-m_I*m_p
m_II = cv2.boxFilter(I*I, -1, (r,r))
var_I = m_II-m_I*m_I
a = cov_Ip/(var_I+eps)
b = m_p-a*m_I
m_a = cv2.boxFilter(a, -1, (r,r))
m_b = cv2.boxFilter(b, -1, (r,r))
return m_a*I+m_b
def getV1(m, r, eps, w, maxV1): #输入rgb图像，值范围[0,1]
'''''计算大气遮罩图像V1和光照值A, V1 = 1-t/A'''
V1 = np.min(m,2) #得到暗通道图像
V1 = guidedfilter(V1, zmMinFilterGray(V1,7), r, eps) #使用引导滤波优化
bins = 2000
ht = np.histogram(V1, bins) #计算大气光照A
d = np.cumsum(ht[0])/float(V1.size)
for lmax in range(bins-1, 0, -1):
if d[lmax]<=0.999:
break
A = np.mean(m,2)[V1>=ht[1][lmax]].max()
V1 = np.minimum(V1*w, maxV1) #对值范围进行限制
return V1,A
def deHaze(m, r=81, eps=0.001, w=0.95, maxV1=0.80, bGamma=False):
Y = np.zeros(m.shape)
V1,A = getV1(m, r, eps, w, maxV1) #得到遮罩图像和大气光照
for k in range(3):
Y[:,:,k] = (m[:,:,k]-V1)/(1-V1/A) #颜色校正
Y = np.clip(Y, 0, 1)
if bGamma:
Y = Y**(np.log(0.5)/np.log(Y.mean())) #gamma校正,默认不进行该操作
return Y
if __name__ == '__main__':
m = deHaze(cv2.imread('test.jpg')/255.0)*255
cv2.imwrite('defog.jpg', m)
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)
return dark
def AtmLight(im,dark):
[h,w] = im.shape[:2]
imsz = h*w
numpx = int(max(math.floor(imsz/1000),1))
darkvec = dark.reshape(imsz)
imvec = im.reshape(imsz,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
def Guidedfilter(im,p,r,eps):
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
mean_II = cv2.boxFilter(im*im,cv2.CV_64F,(r,r))
var_I = mean_II - mean_I*mean_I
a = cov_Ip/(var_I + eps)
b = mean_p - a*mean_I
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)
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__':
import time
fn = './image/15.png'
start = time.time()
for i in range(100):
src = cv2.imread(fn)
I = src.astype('float64')/255
dark = DarkChannel(I,7)
A = AtmLight(I,dark)
te = TransmissionEstimate(I,A,15)
t = TransmissionRefine(src,te)
J = Recover(I,t,A,0.1)
print((time.time()-start)/100)
# cv2.imshow("dark",dark)
# cv2.imshow("t",t)
# cv2.imshow('I',src)
cv2.imshow('J',J)
# cv2.imwrite("./image/J.png",J*255)
cv2.waitKey()
最小值滤波我给用腐蚀来替代了，其实腐蚀就是最小值滤波，最大值滤波是膨胀。这个测试效果还不错。
这份python代码中使用的是暗通道和RGB图像的最小值图像（实际上是一种灰度图）来进行导向滤波，我试着用灰度图和暗通道来做，也是可以的，效果区别不大。
这个python版本跑的还是挺慢的，600-500的图像需要花费近0.1s的时间，c++版本的速度立马提高一倍，代码比python要长一些，就不在这里贴了，相同的图像速度可以提高一倍以上，如果加上GPU加速的话应该可以实现实时处理。
3. 各参数的影响。
暗通道最小值滤波半径r。
这个半径对于去雾效果是有影响的。一定情况下，半径越大去雾的效果越不明显，建议的范围一般是5-25之间，一般选择5,7,9等就会取得不错的效果。w的影响自然也是很大的。
这个值是我们设置的保留雾的程度（c++代码中w是去除雾的程度，一般设置为0.95就可以了）。这个基本不用修改。导向滤波中均值滤波半径。
这个半径建议取值不小于求暗通道时最小值滤波半径的4倍。因为前面最小值后暗通道时一块一块的，为了使得透射率图更加精细，这个r不能过小（很容易理解，如果这个r和和最小值滤波的一样的话，那么在进行滤波的时候包含的块信息就很少，还是容易出现一块一块的形状）。eps，这个值只是保证除号下面不是0，一般取较小，0.001是一个常用的值。
4. notes。
这个去雾算法只针对彩色图像，而且对于低对比度的天空或者水面背景的去雾效果会产生块效应，去雾效果不好，而且这种效应并不能通过调参来避免。
暗通道去雾使得图像整体的亮度有所降低，所以在最后可以自适应的提高亮度来减轻这种现象。
导向滤波在matlab中有现成函数，在opencv contrib里也有函数可以调用，另外为了加速运算可以下采样之后进行滤波然后再上采样恢复。
学习链接：https://www.jianshu.com/p/df9c963a392a