论文原文

清华大学写的文章：基于单幅图像的快速去雾算法，作者：刘倩，陈茂华，周东华。

算法过程

实际上有了这个算法流程就可以写出代码了，不过为了加深理解可以看下面的一些推导。

一些推导

我们知道去雾的步骤主要就是估计全局大气光值A和透射率t(x)，具体可先看这篇文章：https://blog.csdn.net/just_sort/article/details/89634180 。因此，本文就是根据输入图像估计A和L(x)，然后根据雾天退化模型求取去雾后的图像。

估计透射率t(x)

根据上面链接的文章可以知道 $t\left(x\right)\>=1-\frac{H\left(x\right)}{A}$t(x)>=1−AH(x)​…(2)，且论文使用了 $L\left(x\right)$L(x)来代替 $A(1-t(x\left)\right)$A(1−t(x))，即： $L\left(x\right)=A(1-t(x\left)\right)$L(x)=A(1−t(x))…(1)。我们取H(x)三个通道的最小值并记为： $M\left(x\right)=mi{n}_{c\in r,g,b}\left(H^c\right(x\left)\right)$M(x)=minc∈r,g,b​(Hc(x)) …(3)，所以公式2变换为 $t\left(x\right)\>=1-\frac{M\left(x\right)}{A}$t(x)>=1−AM(x)​…(4)，对公式(4)右边进行均值滤波 $meadia{n}_{s_a}(1-\frac{M\left(x\right)}{A})=1-\frac{media{n}_{s_a}\left(M\right(x\left)\right)}{A}=1-\frac{{\sum }_{y\in \Omega \left(x\right)M\left(y\right)}}{A{s}_a^2}$meadiansa​​(1−AM(x)​)=1−Amediansa​​(M(x))​=1−Asa2​∑y∈Ω(x)M(y)​​…(5)。其中 $s_a$sa​代表均值滤波的窗口大小， $\Omega \left(x\right)$Ω(x)表示像素x的 $s_a\times s_a$sa​×sa​的邻域。均值滤波后的结果可以反映 $t\left(x\right)$t(x)的大致趋势，但与真实的 $t\left(x\right)$t(x)还差一定的绝对值，于是，得出透射率的粗略估计值 $t\left(x\right)=1-\frac{M_{ave}\left(x\right)}{A}+\phi \frac{M_{ave}\left(x\right)}{A}=1-\delta \frac{M_{ave}\left(x\right)}{A}$t(x)=1−AMave​(x)​+φAMave​(x)​=1−δAMave​(x)​…(6)。其中 $M_{ave}\left(x\right)=media{n}_{s_a}\left(M\right(x\left)\right),\delta =1-\phi ,\phi \in [0,1]$Mave​(x)=mediansa​​(M(x)),δ=1−φ,φ∈[0,1]，因此 $\phi \in [0,1]$φ∈[0,1]。为了防止去雾后图像出现整体画面偏暗，这里根据图像的均值来调整 $\delta$δ，即 $\delta =\rho m_{ave}$δ=ρmave​,其中 $m_{ave}$mave​是M(x)中所有元素的均值， $\rho$ρ是调节因子。因此可以得到透射率的计算公式 $t\left(x\right)=max(1-min(\rho m_{av},0.9)\frac{M_{ave}\left(x\right)}{A},1-\frac{M\left(x\right)}{A})$t(x)=max(1−min(ρmav​,0.9)AMave​(x)​,1−AM(x)​)…(7)。结合公式(1)推出： $L\left(x\right)=min\left(min\right(\rho m_{av},0.9\left)M_{ave}\right(x),M(x\left)\right)$L(x)=min(min(ρmav​,0.9)Mave​(x),M(x))。

估计全球大气光值

公式5中第一个等式的左侧的表达式取值范围为[0, 1]，由此得出 $A\&gt;=max\left(M_{ave}\right(x\left)\right)$A>=max(Mave​(x))，一般情况下又存在 $A\&lt;=max\left(ma{x}_{c\in r,g,b\left(H^c\right(x\left)\right)}\right)$A<=max(maxc∈r,g,b(Hc(x))​)(KaiMing He的暗通道先验理论)。这样就初步确定了全局大气光的范围，为了能快速获取全局大气光，文章直接取两者的平均值作为全局大气光值，即是：
$A=1/2\left(max\right(H\left(x\right))+max(M_{ave}\left(x\right)\left)\right)$A=1/2(max(H(x))+max(Mave​(x)))…(9)。

代码实现

#include <opencv2/opencv.hpp>
#include <iostream>
#include <algorithm>
#include <vector>
using namespace cv;
using namespace std;

int getMax(Mat src) {
	int row = src.rows;
	int col = src.cols;
	int temp = 0;
	for (int i = 0; i < row; i++) {
		for (int j = 0; j < col; j++) {
			temp = max((int)src.at<uchar>(i, j), temp);
		}
		if (temp == 255) return temp;
	}
	return temp;
}

Mat dehaze(Mat src) {
	double eps;
	int row = src.rows;
	int col = src.cols;
	Mat M = Mat::zeros(row, col, CV_8UC1);
	Mat M_max = Mat::zeros(row, col, CV_8UC1);
	Mat M_ave = Mat::zeros(row, col, CV_8UC1);
	Mat L = Mat::zeros(row, col, CV_8UC1);
	Mat dst = Mat::zeros(row, col, CV_8UC3);
	double m_av, A;
	//get M
	double sum = 0;
	for (int i = 0; i < row; i++) {
		for (int j = 0; j < col; j++) {
			uchar r, g, b, temp1, temp2;
			b = src.at<Vec3b>(i, j)[0];
			g = src.at<Vec3b>(i, j)[1];
			r = src.at<Vec3b>(i, j)[2];
			temp1 = min(min(r, g), b);
			temp2 = max(max(r, g), b);
			M.at<uchar>(i, j) = temp1;
			M_max.at<uchar>(i, j) = temp2;
			sum += temp1;
		}
	}
	m_av = sum / (row * col * 255);
	eps = 0.85 / m_av;
	boxFilter(M, M_ave, CV_8UC1, Size(51, 51));
	double delta = min(0.9, eps*m_av);
	for (int i = 0; i < row; i++) {
		for (int j = 0; j < col; j++) {
			L.at<uchar>(i, j) = min((int)(delta * M_ave.at<uchar>(i, j)), (int)M.at<uchar>(i, j));
		}
	}
	A = (getMax(M_max) + getMax(M_ave)) * 0.5;
	for (int i = 0; i < row; i++) {
		for (int j = 0; j < col; j++) {
			int temp = L.at<uchar>(i, j);
			for (int k = 0; k < 3; k++) {
				int val = A * (src.at<Vec3b>(i, j)[k] - temp) / (A - temp);
				if (val > 255) val = 255;
				if (val < 0) val = 0;
				dst.at<Vec3b>(i, j)[k] = val;
			}
		}
	}
	return dst;
}

int main() {
	Mat src = imread("F:\\fog\\1.jpg");
	Mat dst = dehaze(src);
	cv::imshow("origin", src);
	cv::imshow("result", dst);
	cv::imwrite("F:\\fog\\res.jpg", dst);
	waitKey(0);
	return 0;
}

一些效果图

小结

算法里面有2个参数可以自己调节，滤波的半径和 $\rho$ρ。具体调整方式和取值方式可以参见：http://www.cnblogs.com/Imageshop/p/3410279.html 。欢迎关注我的数字图像处理论文复现工程：https://github.com/BBuf/-Image-processing-algorithm