单尺度Retinex算法学习

最近，看了几篇论文关于retinex算法，其中包括了单尺度SSR与多尺度SSR。今天主要学习一下单尺度SSR，介绍原理及代码，具体实验以后补上。

Retinex算法的基本内容：图像由亮度图像和反射图像构成，分别用L(x,y)和R(x,y)表示。

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

单尺度SSR：

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

                    L(x,y) = I(x,y) * G(x,y)   此处的G(x,y)是高斯卷积函数

       将上述公式转换至对数域：

                   log R = log(I / L) = log I - log L= log I - log(I * G)

实验表明，单尺度Retinex算法可以较好地增强图像。能在保持原始图像亮度的同时压缩图像对比度，并增强图像暗处信息。但当图像中有大块灰度相似的区域时，增强后的图像会产生晕环现象。而这一现象我们可以通过双边滤波方法解决。

代码来自于matlab论坛：http://www.ilovematlab.cn/thread-263388-1-1.html。

clc;clear all ;
I = imread('G:\研二所有论文相关\新海雾图\2013.8.4松花江\松花江\复件 松花江缩放\23.jpg');
figure;imshow(I,[]);title('原图');
I=double(I);
I1=I(:,:,1);
I2=I(:,:,2);
I3=I(:,:,3);
f=double(rgb2gray(I));
%%%%%%%%%%%%%构造高斯核%%%%%%%%%%%%%%%%%%%%%%%
k1=8;
k2=8;
r=161;
% alf=1600;
% alf=161;
alf=20;
nn=floor((r+1)/2);  % nn=81为尺度C，取80最合适
for i=1:r   %从1循环到81
    for j=1:r
%         h(i,j) =exp(-((i-nn)^2+(j-nn)^2)/(k1*alf))/(k2*pi*alf*10000);        %高斯函数
          h(i,j) =exp(-((i-nn)^2+(j-nn)^2)/(2*alf^2)) ;       %高斯函数
   end
end
 h = h/sum(sum(h));%归一化
 %% 把高斯滤波函数转换到频率域
 fgauss=fft2(h,size(h,1),size(h,2));                 %对高斯滤波函数进行傅里叶变换
 fgauss=fftshift(fgauss);                                  %将频域中心移到零点 
 fgauss = ifft2( fgauss);                                  %有变换为空域
%% 利用卷积函数来实现G(x,y)与S(x,y)的卷积效果
q1=conv2(I1,fgauss,'same');  %R通道空域滤波 这里实现的就是Si(x,y)卷积Gx
q2=conv2(I2,fgauss,'same'); %G通道空域滤波
q3=conv2(I3,fgauss,'same');  %B通道空域滤波
q=conv2(f,fgauss,'same');  %灰度
%% 不做对比度拉伸时候的结果 这个结果偏向于蓝色分量偏大的彩色图 得到反射分量log(R(x,y))
G1 = log(I1+1)-log(q1+1);
G2= log(I2+1)-log(q2+1);
G3 = log(I3+1)-log(q3+1);
LLLL=log(f+1)-log(q+1);
GG1=exp(G1);
GG2=exp(G2);
GG3=exp(G3);
ssr=cat(3,G1,G2,G3);
ssr1=cat(3,GG1,GG2,GG3);
figure,imshow(ssr);
figure,imshow(uint8(ssr));
figure,imshow(ssr1);
figure,imshow(uint8(ssr1));
% figure,imshow(mat2gray(uint8(ssr1)));
figure;imshow(LLLL,[]);

ssr=cat(3,abs(G1),abs(G2),abs(G3));  %由于有复数因此这里加abs表示只取实数部分
figure;imshow(mat2gray(ssr));      %显示出了增强效果
figure;imshow(uint8(mat2gray(ssr))); %仍然显示全黑
%%
GG1=exp(G1);
GG2=exp(G2);
GG3=exp(G3);
ssr1=cat(3,abs(GG1),abs(GG2),abs(GG3));
figure;imshow(mat2gray(ssr1));   %这里也显示了增强效果 但是与前面的增强效果差不多
figure;imshow(uint8(ssr1));      %效果很差 很一般