一种图像动态范围调整的自适应方法

注意：这篇文章上次更新于1427天前，文章内容可能已经过时。

~~这也许是我复现的所有论文中效果最好的了吧~~

这也许是我见过最好用的低照度图片了吧🦄

开始之前，背景音乐先安排上。

这篇博客介绍一下论文中的方法，目的不是让除了我以外的人能看懂，因为我也做不到让别人能读懂😓，大概就算是是阅读笔记吧。欲知方法细节，查看论文原文.

方法简介

流程图

算法共分为以下几个步骤：

提取 Y 通道。
根据 Y 通道图的直方图分布情况分别求出两个 $\gamma$ ，使其分别有效作用于暗区和亮区域。
使用 DoG 函数对伽马校正结果图进行局部对比度增强。
对增强结果进行自适应的线性加权融合。
利用 Y 通道的融合结果对图像进行色彩校正得到最终的结果图。

提取 Y 通道

$\begin{aligned} L_{\text {in }}(x, y)=0.299 \cdot I_{\text {in }}^{R}(x, y)+0.587 & \cdot I_{\text {in }}^{G}(x, y) + 0.114 \cdot I_{\text {in }}^{B}(x, y) \end{aligned}$

对 Y 通道进行对数归一化：

$L_{\log }(x, y)=\log _{\alpha}\left(L_{i n}(x, y)+1\right)$

其中：

$\alpha = max(L_{in}(x,y)) + 1$

自适应伽马校正

首先，以 0.5 为阈值将对数归一化后的 Y 通道分为亮区域和暗区域两个部分，即大于 0.5 的像素值被划分在亮区域，小于或等于 0.5 的像素被划分在暗区域。

分别根据亮区域和暗区域的标准差来计算 期望中位数，再通过全局搜索的方式寻找两个 $\gamma$ 使得经过伽马校正后的图像最接近期望中位数。

以下为暗区和亮区期望中位数的计算方法：

$\begin{array}{l} M_{L}=\frac{1}{3}+\sigma_{L} \\ M_{H}=1-\sigma_{H} \end{array}$

其中：L 下标表示暗区；H 下标表示亮区。

两个 $\gamma$ 的全局搜索方法如下：

$\begin{array}{r} \gamma_{L}=\underset{\gamma_{i}}{\arg \min }\left(\left|\operatorname{med}\left(S_{\text {dark }}^{\gamma_{i}}\right)-M_{L}\right|\right) \\ \gamma_{i} \in\{0.1,0.11, \ldots, 1\} \\ \gamma_{H}=\underset{\gamma_{i}}{\arg \min \left(\left|\operatorname{med}\left(S_{\text {bright }}^{\gamma_{i}}\right)-M_{H}\right|\right)} \\ \gamma_{i} \in\{1,1.1, \ldots, 10\} \end{array}$

其中： $med(·)$ 是中位数算子。 $S_{dark}^{\gamma_i}$ 表示的是 $L_{log}(x,y)^{\gamma _i} \le 0.5$ 的像素值。

分别使用以上两个 $\gamma$ 对输入图像进行伽马校正，得到亮区结果和暗区结果。

$\begin{array}{l} L_{d}(x, y)=\left(L_{\log }(x, y)\right)^{\gamma_{L}} \\ L_{b}(x, y)=\left(L_{\log }(x, y)\right)^{\gamma_{H}} \end{array}$

局部对比度增强

使用 DoG 模型分别作用于亮区校正结果图和暗区校正结果图，实现图像的对比度增强。具体方法为：

$\begin{aligned} L_{d}^{\prime}(x, y)=& L_{d}(x, y)+\left(L_{d} * D o G_{\sigma_{c}, \sigma_{s}}\right)(x, y) \\ L_{b}^{\prime}(x, y)=& L_{b}(x, y)+\left(L_{b} * D o G_{\sigma_{c}, \sigma_{s}}\right)(x, y) \\ D o G_{\sigma_{c}, \sigma_{s}}(x, y)=& \frac{1}{2 \pi \sigma_{c}^{2}} \exp \left(-\frac{x^{2}+y^{2}}{2 \sigma_{c}^{2}}\right) -\frac{1}{2 \pi \sigma_{s}^{2}} \exp \left(-\frac{x^{2}+y^{2}}{2 \sigma_{s}^{2}}\right) \end{aligned}$

自适应融合

这里使用一种依赖于亮度的自适应线性加权融合方法。

其权重的计算方式为：

$w(x, y)=\exp \left(-\frac{L_{b}(x, y)^{2}}{2 \sigma_{w}^{2}}\right)$

此处设置 $\sigma_w = 0.5$ 。

融合方法为：

$L_{\text {out }}(x, y)=w(x, y) \cdot L_{d}^{\prime}(x, y)+(1-w(x, y)) \cdot L_{b}^{\prime}(x, y)$

另外，此处对输出结果进行了一步线性归一化操作。

自适应的色彩恢复

此处依然是一种亮度依赖的自适应方法。

$\begin{aligned} I_{\text {out }}^{c}(x, y) &=L_{\text {out }}(x, y)\left(\frac{I_{i n}^{c}(x, y)}{L_{\text {in }}(x, y)}\right)^{s(x, y)} \\ s(x, y) &=1-\tanh \left(L_{b}(x, y)\right) \end{aligned}$

Matlab实现

分区

分区操作只需要记录亮区像素和暗区像素的值就可以，因此此处采用遍历图像的方式对每个像素点进行判断，大于0.5 则加入亮区数组，否则加入暗区数组。

function [Ld,Lb] = fen_qu(Lin)
i = 1;
j = 1;
[m,n] = size(Lin);
Ld = 0;
Lb = 0;
for k = 1:m*n
    x = Lin(k);
    if x <= 0.5
        Ld(i) = x;
        i = i + 1;
    end
    if x > 0.5
        Lb(j) = x;
        j = j + 1;
    end
end
end

DoG 算子

DoG 算子实际上就是两个高斯函数的差，因此这里先实现一个二维的高斯函数。

二维高斯函数也可以使用 Matlab 内置函数fspecial来实现。

function gaussFun = GaussFun2D(hsize) 

% sigma = double(hsize)/6.0; % 3*sigma = radius, i.e., 6*sigma = diameter (hsize)
% gaussFun = fspecial('gaussian',hsize,sigma);

gaussFun = zeros(hsize,hsize);
rr = (hsize+1)/2;
%rr = (hsize-1)/2;
for i=1:hsize
    for j=1:hsize
				tmpV = (i-rr)*(i-rr) + (j-rr)*(j-rr);
				gaussFun(i,j) = exp(-9.0*tmpV/rr/rr);
				%gaussFun(i,j) = exp(-4.5*tmpV/rr/rr);
				%sigma = rr/3;
				%gaussFun(i,j) = exp(-tmpV/(2*sigma*sigma))/(2*pi*sigma*sigma);
    end
end
end

DoG 算子的实现如下：

function outImg = DoG(Img,rr_RF,rr_IF,m_qA1,m_qA2,antImg)

th1 = 2*rr_RF + 1; 
hsize = th1;
pa_gauss1 = GaussFun2D(hsize);
for i=1:th1
    for j=1:th1
        di = i-rr_RF-1;
        dj = j-rr_RF-1;
        tmpV = di*di + dj*dj;
        if tmpV>rr_RF*rr_RF
            pa_gauss1(i,j) = 0.0;
        end
    end
end
pa_gauss1 = m_qA1*pa_gauss1/sum(sum(pa_gauss1));


th2  = 2*rr_IF + 1;
hsize = th2;
pa_gauss2 = GaussFun2D(hsize);
for i=1:th2
    for j=1:th2
        di = i-rr_IF-1;
        dj = j-rr_IF-1;
        tmpV = di*di + dj*dj;
        if tmpV>rr_IF*rr_IF || tmpV<=rr_RF*rr_RF
            pa_gauss2(i,j) = 0.0;
        end
    end
end
pa_gauss2 = m_qA2*pa_gauss2/sum(sum(pa_gauss2));

SurroundImgLists = imfilter(antImg,pa_gauss2,'conv','symmetric');
CenterImgLists = imfilter(Img,pa_gauss1,'conv','symmetric');

outImg = CenterImgLists - SurroundImgLists;
outImg(outImg<0.0) = 0.0;
end

线性归一化

function [outputs] = Normalization(arg)
    MAX = max(arg(:));
    MIN = min(arg(:));
    outputs = (arg - MIN)./ (MAX - MIN);
end

其他操作

最后，将其他操作写在主函数中，详见代码注释。

function outputs = main(input)
I = im2double(input);
R = I(:,:,1);
G = I(:,:,2);
B = I(:,:,3);

% 计算 Y 通道
Lin = 0.299*R + 0.587*G +0.114*B;

% 对数归一化
alpha = max(max(Lin)) + 1;
Llog = log(Lin + 1) / log(alpha);

% 计算期望中位数
[Ld , Lb] = fen_qu(Llog);
ML = 1/3 + std(Ld);
MH = 1 - std(Lb);

% 搜索最优 γ
Min1 = zeros(1,91);
Min2 = zeros(1,91);
k = 1;
for i = 0.1:0.01:1
    clear Ld
    [Ld , ~] = fen_qu(Llog.^i);
    Min1(k) = abs(median(Ld) - ML);
    k = k + 1;
end
k = 1;
for i = 1:0.1:10
    clear Lb
    [~ , Lb] = fen_qu(Llog.^i);
    Min2(k) = abs(median(Lb) - MH);
    k = k + 1;
end
gamma_L = (find(Min1 == min(Min1)) - 1) * 0.01 + 0.1;
gamma_H = (find(Min2 == min(Min2)) - 1) * 0.1 + 1;

% 伽马校正
Ld = Llog .^ gamma_L;
Lb = Llog .^ gamma_H;

% 局部对比度增强
Ld_p = Ld + DoG(Ld,1,3,1,1,Ld);
Lb_p = Lb + DoG(Lb,1,3,1,1,Lb);

% 自适应权重计算与图像融合
w = exp(-(Lb.^2)./(2*0.5*0.5));
L_out = w .* Ld_p + (1-w) .* Lb_p;
L_out = Normalization(L_out);

% 自适应色彩恢复
s = 1 - tanh(Lb);

R_out = L_out .* ((R ./ Lin).^s);
G_out = L_out .* ((G ./ Lin).^s);
B_out = L_out .* ((B ./ Lin).^s);

% 输出
outputs = cat(3,R_out,G_out,B_out);
end

几张低照度图像

左为原图，右为结果图。

703

785

几张水下图像

左为原图，右为结果图。

184

910