Retinex学习笔记（一）

Retinex学习

综述

理论

Retinex图像增强的基本思想：去除照射光影响，保留物体自身的反射属性。

核心思想：保留图像细节信息的前提下，调整图像的对比度和亮度。

Retinex 理论认为图像 I(x, y)是由照度图像与反射图像组成。前者指的是物体的入射分量的信息，用 L(x, y) 表示；后者指的是物体的反射部分，用 R(x, y) 表示。公式：

I(x, y) =R(x, y) * L(x, y)

同时，由于对数形式与人类在感受亮度的过程属性最相近，因此将上述过程转换到对数域进行处理，这样做也将复杂的乘法转换为加法： i(x, y) = r(x, y) + l(x, y)

照度图和反射图

照度图和反射图的区别是，照度图表示光源对某个表面的照射强度，而反射图表示物体表面对光线的反射特性。照度图和反射图都可以用于计算光照和反射效果，但是它们的输入和输出不同。

照度图的输入是光源的位置、强度和颜色，输出是表面上不同位置的照度值。反射图的输入是物体表面的材质、颜色和法线，输出是表面上不同位置的反射率和颜色

流程

高斯滤波：对输入图像进行高斯模糊处理，以消除图像中的噪声和细节**。

计算照度图像：将处理后的图像转换为对数空间，并通过高斯滤波计算出照度图像，即物体的入射分量信息。

计算反射图像：通过计算原始图像和照度图像的差异得到反射图像，即物体的反射部分。

反对数变换：将照度图像和反射图像转换回线性空间。

图像拉伸：对输出图像进行亮度调整，以增强图像的对比度和细节。

发展历程

Frankle-McCann Retinex

原理理解

Frankle-McCann Retinex是一种用于图像增强的算法，它可以分离图像中的光照和反射成分，从而改善图像的对比度和色彩。它是基于Retinex理论的，该理论认为人眼可以根据场景中的光线变化来调整视觉感知。 Frankle-McCann Retinex算法的基本思想是，从图像中选取一个参考点，然后递归地计算与参考点相邻的像素点的亮度和颜色，直到覆盖整个图像。这个过程可以重复多次，以达到更好的效果。

主要就是利用一个对数公式计算当前像素点的亮度，这个亮度与周围像素点的亮度有关。

代码理解

function Test()
    % 读取图像
    ImOriginal = imread('fig5.tif');
​
    % 获取图像尺寸
    [m, n, z] = size(ImOriginal);
​
    % 初始化输出图像
    ImOut = zeros(m, n, z);
​
    % 对每个通道进行Retinex增强
    for i = 1:z
        ImChannel = log(double(ImOriginal(:, :, i)) + eps);
        ImOut(:, :, i) = retinex_frankle_mccann(ImChannel, 4); % 使用Frankle-McCann Retinex算法增强
        ImOut(:, :, i) = exp(ImOut(:, :, i));
        a = min(min(ImOut(:, :, i)));
        b = max(max(ImOut(:, :, i)));
        ImOut(:, :, i) = ((ImOut(:, :, i) - a) / (b - a)) * 255; % 将增强后的像素值映射到0~255之间
    end
​
    % 转换图像数据类型为uint8
    ImOut = uint8(ImOut);
​
    % 显示原始图像和增强后的图像
    figure(1);
    imshow(ImOriginal);
    title('原始图像');
​
    figure(2);
    imshow(ImOut);
    title('Retinex增强后的图像');
​
    % 将增强后的图像保存为tt.tif
    imwrite(ImOut, 'tt.tif');
​
end
​
function Retinex = retinex_frankle_mccann(L, nIterations)
    global RR IP OP NP Maximum
    RR = L;
    Maximum = max(L(:)); % 获取图像中的最大像素值
    [nrows, ncols] = size(L);
​
    shift = 2^(fix(log2(min(nrows, ncols))) - 1); % 初始位移
    OP = Maximum * ones(nrows, ncols); % 初始化Old Product
​
    % 迭代Retinex增强
    while (abs(shift) >= 1)
        for i = 1:nIterations
            CompareWith(0, shift); % 水平方向的比较
            CompareWith(shift, 0); % 垂直方向的比较
        end
        shift = -shift / 2; % 更新位移
    end
​
    Retinex = NP; % 返回增强后的图像
​
end
​
function CompareWith(s_row, s_col)
    global RR IP OP NP Maximum
    IP = OP;
    if (s_row + s_col > 0)
        % 比较并更新像素值
        IP((s_row + 1):end, (s_col + 1):end) = OP(1:(end - s_row), 1:(end - s_col)) + ...
            RR((s_row + 1):end, (s_col + 1):end) - RR(1:(end - s_row), 1:(end - s_col));
    else
        % 比较并更新像素值
        IP(1:(end + s_row), 1:(end + s_col)) = OP((1 - s_row):end, (1 - s_col):end) + ...
            RR(1:(end + s_row), 1:(end + s_col)) - RR((1 - s_row):end, (1 - s_col):end);
    end
​
    IP(IP > Maximum) = Maximum; % 重置像素值，使其不超过最大值
    NP = (IP + OP) / 2; % 取新旧像素值的平均值作为更新后的像素值
    OP = NP; % 准备下一轮迭代比较
​
end
​

1. Test函数：这是主函数，用于读取图像、执行Retinex增强，并显示和保存结果。

2. retinex_frankle_mccann函数：这是Retinex算法的核心函数，对输入图像进行增强。

3. CompareWith函数：辅助函数，用于在Retinex算法的迭代过程中进行图像的水平和垂直比较。

参数调整

想调节参数以获得更好的Retinex增强效果，主要可以考虑调整以下部分：

1. retinex_frankle_mccann函数中的nIterations参数： nIterations控制Retinex增强算法的迭代次数。增加nIterations可以增强图像的对比度和细节，但可能会导致过度增强和图像噪声增加。减少nIterations可能会导致增强效果不够显著。你可以尝试不同的值来找到最佳的增强效果。

2. 在主函数Test中，对于Retinex增强前后的像素值映射：

   scssCopy code
   ImOut(:, :, i) = ((ImOut(:, :, i) - a) / (b - a)) * 255;

   上述代码将像素值映射到0到255之间。你可以根据实际需求，调整映射范围或使用其他映射方法，以便更好地显示增强后的图像。

3. 在CompareWith函数中，有一个用于控制像素值重置的操作：

   scssCopy code
   IP(IP > Maximum) = Maximum;

   如果你觉得像素值重置过于严格，你可以调整Maximum的值，或者尝试其他像素值截断或缩放的方法。

4. 对于输入图像的预处理部分：

   cCopy code
   ImChannel = log(double(ImOriginal(:, :, i)) + eps);

   上述代码使用了对数变换对输入图像进行预处理。你可以尝试使用其他预处理方法，比如伽马校正或直方图均衡化，来改善Retinex增强效果。

McCann99 Retinex

原理

McCann99 Retinex 算法是一种基于Retinex理论的图像增强算法，它利用图像金字塔建立多分辨率描述，自顶向下逐层迭代，提高增强效率。它可以改善图像的对比度和色彩。但是，它对输入图像的长宽有严格的限制，需要满足一定的条件24。有些改进型算法对这个问题进行了优化。

最底层为原始图像，最上面为设定好的下采样（池化）图像：n*m，每下一层大小为上一层的2倍，所以对于原始图像有比较严格的限制。

优点： 颜色恒常性、动态范围压缩大（像素点丰富）、色彩逼真度高（图像高保真）

缺点： 增强后的图像存在“光晕”现象，即在图像色彩交界处渐变 光晕缺陷同样也存在于其他 Retinex算法之中

代码

function Test()
    % 读取图像
    ImOriginal = imread('fig5.tif');
​
    % 获取图像尺寸
    [m, n, z] = size(ImOriginal);
​
    % 初始化输出图像
    ImOut = zeros(m, n, z);
​
    % 对每个通道进行Retinex增强
    for i = 1:z
        % 使用对数变换预处理图像
        ImChannel = log(double(ImOriginal(:,:,i))+eps);
        % 使用Frankle-McCann Retinex算法增强
        ImOut(:,:,i) = retinex_mccann99(ImChannel, 4); 
        % 对增强后的像素值进行映射
        ImOut(:,:,i) = exp(ImOut(:,:,i));
        a = min(min(ImOut(:,:,i)));
        b = max(max(ImOut(:,:,i)));
        ImOut(:,:,i) = ((ImOut(:,:,i) - a) / (b - a)) * 255;     
    end
​
    % 转换图像数据类型为uint8
    ImOut = uint8(ImOut);
​
    % 显示原始图像和增强后的图像
    figure(1);
    imshow(ImOriginal);
    title('原始图像');
​
    figure(2);
    imshow(ImOut);
    title('Retinex增强后的图像');
​
end
​
function Retinex = retinex_mccann99(L, nIterations)
    % 输入: L - 对数变换后的单通道强度图像
    %      nIterations - Retinex算法的迭代次数
​
    global OPE RRE Maximum
    [nrows, ncols] = size(L); % 获取输入图像的尺寸
    nLayers = ComputeLayers(nrows, ncols); % 计算金字塔的层数
    nrows = nrows / (2^nLayers); % 计算第0层的图像大小
    ncols = ncols / (2^nLayers);
    
    if (nrows * ncols > 25) % 图像面积大于25不进行处理
        error('invalid image size.')
    end
    
    Maximum = max(L(:)); % 获取图像中的最大像素值
    OP = Maximum * ones([nrows ncols]); % 初始化Old Product
    
    for layer = 0:nLayers
        RR = ImageDownResolution(L, 2^(nLayers-layer)); % 对输入图像降采样到所需层大小
        OPE = [zeros(nrows, 1) OP zeros(nrows, 1)]; % 使用零进行填充以增加OP的列数
        OPE = [zeros(1, ncols+2); OPE; zeros(1, ncols+2)]; % 使用零进行填充以增加OP的行数
        RRE = [RR(:, 1) RR RR(:, end)]; % 使用RR的列值进行填充以增加RRE的列数
        RRE = [RRE(1, :); RRE; RRE(end, :)]; % 使用RR的行值进行填充以增加RRE的行数
        
        for iter = 1:nIterations
            CompareWithNeighbor(-1, 0); % 对八个方向的邻域进行比较和更新
            CompareWithNeighbor(-1, 1);
            CompareWithNeighbor(0, 1);
            CompareWithNeighbor(1, 1);
            CompareWithNeighbor(1, 0);
            CompareWithNeighbor(1, -1);
            CompareWithNeighbor(0, -1);
            CompareWithNeighbor(-1, -1);
        end
        
        NP = OPE(2:(end-1), 2:(end-1));
        OP = NP(:, [fix(1:0.5:ncols) ncols]); % 下采样后的像素放大到原始大小
        OP = OP([fix(1:0.5:nrows) nrows], :);
        nrows = 2 * nrows; ncols = 2 * ncols; % 将图像大小扩大两倍，准备处理下一层
        
    end
    
    Retinex = NP; % 返回增强后的图像
    
end
​
function CompareWithNeighbor(dif_row, dif_col)
    % 输入: dif_row, dif_col - 八个方向的差值
    global OPE RRE Maximum
    
    % Ratio-Product运算
    IP = OPE(2 + dif_row:(end - 1 + dif_row), 2 + dif_col:(end - 1 + dif_col)) + ...
         RRE(2:(end - 1), 2:(end - 1)) - RRE(2 + dif_row:(end - 1 + dif_row), 2 + dif_col:(end - 1 + dif_col));
    IP(IP > Maximum) = Maximum; % 超过最大值的像素值重置为最大值
    
    % 忽略邻域不可定义的行或列得到的结果
    if (dif_col == -1) IP(:, 1) = OPE(2:(end - 1), 2); end
    if (dif_col == +1) IP(:, end) = OPE(2:(end - 1), end - 1); end
    if (dif_row == -1) IP(1, :) = OPE(2, 2:(end - 1)); end
    if (dif_row == +1) IP(end, :) = OPE(end - 1, 2:(end - 1)); end
    
    NP = (OPE(2:(end - 1), 2:(end - 1)) + IP) / 2; % 平均操作
    OPE(2:(end - 1), 2:(end - 1)) = NP; % 更新OPE
    
end
​
function Layers = ComputeLayers(nrows, ncols)
    % 计算层数，最大公约数为2的整数次幂
    power = 2^fix(log2(gcd(nrows, ncols)));
    while(power > 1 && ((rem(nrows, power) ~= 0) || (rem(ncols, power) ~= 0)))
        power = power/2;
    end
    Layers = log2(power);
    
end
​
function Result = ImageDownResolution(A, blocksize)
    % 输入: A - 输入矩阵
    %      blocksize - 下采样的块大小
    % 输出: Result - 下采样后的结果矩阵
    
    [rows, cols] = size(A);
    result_rows = rows / blocksize;
    result_cols = cols / blocksize;
    Result = zeros([result_rows result_cols]);
    
    for crt_row = 1:result_rows
        for crt_col = 1:result_cols
            % 每个像素点通过块的均值得到
            Result(crt_row, crt_col) = mean2(A(1+(crt_row-1)*blocksize:crt_row*blocksize, ...
                                               1+(crt_col-1)*blocksize:crt_col*blocksize));
        end
    end
    
end
​