MSR MATLAB程序

clc
% clear all
I=imread（‘281.bmp‘）;
I=im2double（I）;
%将灰度值的数据类型由BYTE转换为Double型 
[m n bit]=size（I）;

R=fftshift（fft2（I（::1）））;
G=fftshift（fft2（I（::2）））;
B=fftshift（fft2（I（::3）））;

s=3;
sigma=[15 75 115];
for i=1:s
 H（::i）=fftshift（fft2（hfunction（mnsigma（i））））;
end

D=zeros（mn）;
for i=1:m
    for j=1:n
        D（ij）=（-1）^（i+j）;
    end
end

for i=1:s
    r（::i）=abs（ifft2（ifftshift（R.*H（::i）.*D）））;
    g（::i）=abs（ifft2（ifftshift（G.*H（::i）.*D）））;
    b（::i）=abs（ifft2（ifftshift（B.*H（::i）.*D）））;
end

w=[1.0 1.0 1.0]/s;
rr=zeros（mn）;
gg=zeros（mn）;
bb=zeros（mn）;
for i=1:s
    rr=rr+w（i）.*（log2（I（::1）+0.0001）-log2（r（::i）））;
    gg=gg+w（i）.*（log2（I（::2）+0.0001）-log2（g（::i）））;
    bb=bb+w（i）.*（log2（I（::3）+0.0001）-log2（b（::i）））;
end

I1=zeros（size（I））;
I1（::1）=rr/max（rr（:））;
I1（::2）=gg/max（gg（:））;
I1（::3）=bb/max（bb（:））;

avg=mean（mean（mean（I1）））;
div=0;
for i=1:m
    for j=1:n
        for q=1:bit
            div=div+（I1（ijq）-avg）^2;
        end
    end
end
div=sqrt（div/（m*n*bit））;

k=3;
fmax=avg+k*div;
%表示增强后图像中所有像素点灰度值中的最大值
fmin=avg-k*div;
%表示增强后图像中所有像素点灰度值中的最小值 
for i=1:m
    for j=1:n
        for q=1:bit
            if I1（ijq）                I1（ijq）=fmin;
            end
            if I1（ijq）>fmax
                I1（ijq）=fmax;
            end
        end
    end
end
idensity=mean（mean（mean（I）））;
I1=（I1+（idensity/2）-fmin）/（fmax+（idensity/2）-fmin）;
%对像素点灰度值校正后的图象进行线形拉伸
figure
imshow（I）;
title（‘原始图像‘）
figure
imshow（I1）;

B=mean（I3）;
nz=find（B）;
B（nz） =1./B（nz）;
%linear
I2=zeros（size（I））;
I2（::1）=rr.*（I（::1）.*B）;
I2（::2）=gg.*（I（::2）.*B）;
I2（::3）=bb.*（I（::3）.*B）;
I2（::1）=I2（::1）/max（max（max（I2（::1））））;
I2（::2）=I2（::2）/max（max（max（I2（::2））））;
I2（::3）=I2（::3）/max（max（max（I2（::3））））;

avg2=mean（mean（mean（I2）））;
div2=0;
for i=1:m
    for j=1:n
        for q=1:bit
            div2=div2+（I2（ijq）-avg2）^2;
        end
    end
end
div2=sqrt（div2/（m*n*bit））;

k2=3;
fmax2=avg2+k2*div2;
fmin2=avg2-k2*div2;
for i=1:m
    for j=1:n
        for q=1:bit
            if I2（ijq）                I2（ijq）=fmin2;
            end
            if I2（ijq）>fmax2
                I2（ijq）=fmax2;
            end
        end
    end
end
I2=（I2+（idensity/2）-fmin2）/（fmax2+（idensity/2）-fmin2）;
figureimshow（I2）;
%log 
I3=zeros（size（I））;
I3（::1）=rr.*log2（1+I（::1）.*B）;
I3（::2）=gg.*log2（1+I（::2）.*B）;
I3（::3）=bb.*log2（1+I（::3）.*B）;
I3（::1）=I3（::1）/max（max（max（I3（::1））