最近,看了几篇论文关于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)); %效果很差 很一般