论文原文
清华大学写的文章:基于单幅图像的快速去雾算法,作者:刘倩,陈茂华,周东华。
算法过程
实际上有了这个算法流程就可以写出代码了,不过为了加深理解可以看下面的一些推导。
一些推导
我们知道去雾的步骤主要就是估计全局大气光值A和透射率t(x),具体可先看这篇文章:https://blog.csdn.net/just_sort/article/details/89634180 。因此,本文就是根据输入图像估计A和L(x),然后根据雾天退化模型求取去雾后的图像。
估计透射率t(x)
根据上面链接的文章可以知道 t(x)>=1−AH(x)…(2),且论文使用了 L(x)来代替 A(1−t(x)),即: L(x)=A(1−t(x))…(1)。我们取H(x)三个通道的最小值并记为: M(x)=minc∈r,g,b(Hc(x)) …(3),所以公式2变换为 t(x)>=1−AM(x)…(4),对公式(4)右边进行均值滤波 meadiansa(1−AM(x))=1−Amediansa(M(x))=1−Asa2∑y∈Ω(x)M(y)…(5)。其中 sa代表均值滤波的窗口大小, Ω(x)表示像素x的 sa×sa的邻域。均值滤波后的结果可以反映 t(x)的大致趋势,但与真实的 t(x)还差一定的绝对值,于是,得出透射率的粗略估计值 t(x)=1−AMave(x)+φAMave(x)=1−δAMave(x)…(6)。其中 Mave(x)=mediansa(M(x)),δ=1−φ,φ∈[0,1],因此 φ∈[0,1]。为了防止去雾后图像出现整体画面偏暗,这里根据图像的均值来调整 δ,即 δ=ρmave,其中 mave是M(x)中所有元素的均值, ρ是调节因子。因此可以得到透射率的计算公式 t(x)=max(1−min(ρmav,0.9)AMave(x),1−AM(x))…(7)。结合公式(1)推出: L(x)=min(min(ρmav,0.9)Mave(x),M(x))。
估计全球大气光值
公式5中第一个等式的左侧的表达式取值范围为[0, 1],由此得出 A>=max(Mave(x)),一般情况下又存在 A<=max(maxc∈r,g,b(Hc(x)))(KaiMing He的暗通道先验理论)。这样就初步确定了全局大气光的范围,为了能快速获取全局大气光,文章直接取两者的平均值作为全局大气光值,即是:
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个参数可以自己调节,滤波的半径和 ρ。具体调整方式和取值方式可以参见:http://www.cnblogs.com/Imageshop/p/3410279.html 。欢迎关注我的数字图像处理论文复现工程:https://github.com/BBuf/-Image-processing-algorithm