單尺度Retinex（matlab）

說明：此程式碼是針對灰度影像的

clc;clear;
I=imread('DR_ck_0001.tiff');
 
% Ir=I(:,:,1);%提取紅色分量
% Ig=I(:,:,2);%提取綠色分量
% Ib=I(:,:,3);%提取藍色分量
% I1=SSR(Ir); %自定義單尺度Retinex處理函式
% I2=SSR(Ig);
% I3=SSR(Ib);
% In=cat(3,I1,I2,I3);          %cat用於構造多維陣列
In = SSR(I);
t = exp(-10);
enhanced =  log(double(I)+t)-log(In+t);
enhanced = normal(enhanced);
subplot(1,3,1);imshow(I);
title('A). 原始影像');  
subplot(1,3,2);imshow(In,[]);
title('B). 單尺度Retinex');
subplot(1,3,3);imshow(enhanced,[]);
title('C). 單尺度Retinex');
 
%——————————————劃重點！！！！！！！！！！！！！！！！！！！！————————————————
 
function G = SSR(I)
%G=I;
[m,n]=size(I);
I=double(I);
x0=floor((m+1)/2);y0=floor((n+1)/2);
c = 70;
K=1/(sqrt(2*pi)*c);
for x=1:m
    for y=1:n
            G(x,y)=K*(exp(-((x)^2+(y)^2)/(2*c^2))); %高斯函式
    end
end
 
Gfft=fft2(G);
Ifft=fft2(I);
Lfft=Gfft.*Ifft;
L=ifft2(Lfft);
 
r=log(I+1)-log(L+1);                   %加1是為了防止對數為0時，log0沒有定義
 
%線性拉伸，將灰度範圍轉換到0-255，（直接求反對數的效果不好）
MIN=min(min(r));
MAX=max(max(r));
new=((r-MIN)*255/(MAX-MIN));
A=new;
end
function out =normal(img)
maxx = max(max(img));
minn = min(min(img));
out=(img-minn)./(maxx-minn);
end