Retinex image de blur (including MATLAB code)

Retinex image de blur

Retinex theory:
Retinax is composed of Retina + Cortex, which is called the theory of retina cortex. The theory points out that the information that an object can be observed is determined by two factors: the reflection property of the object itself and the light intensity around the object. Among them, the intensity of light determines the dynamic range of all pixels in the original image, and the inherent attribute (color) of the original image is determined by the reflection coefficient of the object itself. Therefore, this method divides an image into two different images: reflection image and brightness image, then removes the influence of light, and retains the inherent attributes of the object.

At present, this method has three main branches: single scale Retinex (SSR), multi-scale Retinex (MSR) and multi-scale Retinex with color recovery (MSRCR).
The basic idea of Retinex theory: in the original image, some methods are used to remove or reduce the impact of the incident image, so as to retain the essential reflection attribute image of the object.

Single scale Retinax (SSR):

Algorithm idea:
1. Read the original image, take out the numerical matrix of RGB three channels, and convert it into double type.
2. The convolution of each channel is processed by Gauss template.
The RGB image and the convoluted image are transformed into the logarithmic domain, and then the RGB image is subtracted, and then the inverse number is transformed into the real number domain.
4. The image of each channel obtained is stretched linearly, and the processed image can be combined with RGB.
Algorithm implementation (MATLAB):

close all; clear all; clc
I = imread('E:\Underwater target recognition\chinamm2019uw_train\jpegimage\000000.jpg');
I_r = double(I(:,:,1));
I_g = double(I(:,:,2));
I_b = double(I(:,:,3));

I_r_log = log(I_r+1);
I_g_log = log(I_g+1);
I_b_log = log(I_b+1);

Rfft1 = fft2(I_r);
Gfft1 = fft2(I_g);
Bfft1 = fft2(I_b);

%  SSR algorithm
[m,n] = size(I_r);
sigma = 200;
f = fspecial('gaussian', [m, n], sigma);
efft1 = fft2(double(f));

D_r = ifft2(Rfft1.*efft1);
D_g = ifft2(Gfft1.*efft1);
D_b = ifft2(Bfft1.*efft1);

D_r_log = log(D_r + 1);
D_g_log = log(D_g + 1);
D_b_log = log(D_b + 1);

R = I_r_log - D_r_log;
G = I_g_log - D_g_log;
B = I_b_log - D_b_log;

R = exp(R);
MIN = min(min(R)); 
MAX = max(max(R));
R = (R - MIN)/(MAX - MIN);
R = adapthisteq(R);
G = exp(G);
MIN = min(min(G)); 
MAX = max(max(G));
G = (G - MIN)/(MAX - MIN);
G = adapthisteq(G);
B = exp(B);
MIN = min(min(B)); 
MAX = max(max(B));
B = (B - MIN)/(MAX - MIN);
B = adapthisteq(B);

J = cat(3, R, G, B);

figure;
subplot(121);imshow(I);
subplot(122);imshow(J);
figure;imshow(J)

Result:

Multiscale Retinex (MSR):

Because SSR needs to pursue a perfect balance between color fidelity and details, this balance is not suitable. MSR is proposed to solve this problem. It aims at one image to be processed by Gauss filter on different scales, and then the images under multiple scales are weighted.
It is actually a stack operation of SSR.
Code:

close all; clear all; clc
I = imread('E:\Underwater target recognition\chinamm2019uw_train\jpegimage\000000.jpg');
I_r = double(I(:,:,1));
I_g = double(I(:,:,2));
I_b = double(I(:,:,3));

I_r_log = log(I_r+1);
I_g_log = log(I_g+1);
I_b_log = log(I_b+1);

Rfft1 = fft2(I_r);
Gfft1 = fft2(I_g);
Bfft1 = fft2(I_b);[m,n] = size(I_r);
sigma1 = 15;
sigma2 = 80;
sigma3 = 200;
f1 = fspecial('gaussian', [m, n], sigma1);
f2 = fspecial('gaussian', [m, n], sigma2);
f3 = fspecial('gaussian', [m, n], sigma3);
efft1 = fft2(double(f1));
efft2 = fft2(double(f2));
efft3 = fft2(double(f3));

D_r1 = ifft2(Rfft1.*efft1);
D_g1 = ifft2(Gfft1.*efft1);
D_b1 = ifft2(Bfft1.*efft1);
D_r_log1 = log(D_r1 + 1);
D_g_log1 = log(D_g1 + 1);
D_b_log1 = log(D_b1 + 1);
R1 = I_r_log - D_r_log1;
G1 = I_g_log - D_g_log1;
B1 = I_b_log - D_b_log1;

D_r2 = ifft2(Rfft1.*efft2);
D_g2 = ifft2(Gfft1.*efft2);
D_b2 = ifft2(Bfft1.*efft2);
D_r_log2 = log(D_r2 + 1);
D_g_log2 = log(D_g2 + 1);
D_b_log2 = log(D_b2 + 1);
R2 = I_r_log - D_r_log2;
G2 = I_g_log - D_g_log2;
B2 = I_b_log - D_b_log2;

D_r3 = ifft2(Rfft1.*efft3);
D_g3 = ifft2(Gfft1.*efft3);
D_b3 = ifft2(Bfft1.*efft3);
D_r_log3 = log(D_r3 + 1);
D_g_log3 = log(D_g3 + 1);
D_b_log3 = log(D_b3 + 1);
R3 = I_r_log - D_r_log3;
G3 = I_g_log - D_g_log3;
B3 = I_b_log - D_b_log3;

R = 0.1*R1 + 0.4*R2 + 0.5*R3;
G = 0.1*G1 + 0.4*G2 + 0.5*G3;
B = 0.1*B1 + 0.4*B2 + 0.5*B3;

R = exp(R);
MIN = min(min(R)); 
MAX = max(max(R));
R = (R - MIN)/(MAX - MIN);
R = adapthisteq(R);
G = exp(G);
MIN = min(min(G)); 
MAX = max(max(G));
G = (G - MIN)/(MAX - MIN);
G = adapthisteq(G);
B = exp(B);
MIN = min(min(B)); 
MAX = max(max(B));
B = (B - MIN)/(MAX - MIN);
B = adapthisteq(B);

J = cat(3, R, G, B);

figure;
subplot(121);imshow(I);
subplot(122);imshow(J,[]);

figure;imshow(J)

Result:

Multiscale Retinex with color recovery (MSRCR):

Because of the color difference in the process of using SSR and MSR, MSRCR adjusts the result of MSR according to a certain proportion to restore the original proportion value, specifically by introducing the color recovery factor C, to make up for the defect of color distortion caused by the enhancement of contrast in the local area of the image.
Code:

% close all; clear all; clc
% I = imread('E:\Underwater target recognition\chinamm2019uw_train\jpegimage\000000.jpg');
% figure;imshow(I)
% I_r = I(:,:,1);
% I_g = I(:,:,2);
% I_b = I(:,:,3);
% figure;
% subplot(131);imshow(I_r);
% subplot(132);imshow(I_g);
% subplot(133);imshow(I_b);
% 
% I_r = double(I_r);
% I_g = double(I_g);
% I_b = double(I_b);
% 
% [m, n] = size(I_r);
% Rlog = log(I_r+1);
% Rfft2 = fft2(I_r);
% 
% Glog = log(I_g+1);
% Gfft2 = fft2(I_g);
% 
% Blog = log(I_b+1);
% Bfft2 = fft2(I_b);
% 
% sigma1 = 128;
% f = fspecial('gaussian', [m,n], sigma1);
% Efft1 = fft2(double(f));
% 
% DI_r = Rfft2.* Efft1;
% DR = ifft2(DI_r);
% DRlog = log(DR+1);
% Rr1 = Rlog - DRlog;
% figure;imshow(Rr1);
% 
% DI_g = Gfft2.* Efft1;
% DG = ifft2(DI_g);
% DGlog = log(DG+1);
% Gg1 = Glog - DGlog;
% figure;imshow(Gg1);
% 
% DI_b = Bfft2.* Efft1;
% DB = ifft2(DI_b);
% DBlog = log(DB+1);
% Bb1 = Blog - DBlog;
% figure;imshow(Bb1);
% 
% Rr1 = Rr1 .* 255/max(max(Rr1));
% Gg1 = Gg1 .* 255/max(max(Gg1));
% Bb1 = Bb1 .* 255/max(max(Bb1));
% 
% II(:,:,1) = Rr1;
% II(:,:,2) = Gg1;
% II(:,:,3) = Bb1;
% figure;
% subplot(121);imshow(I);
% subplot(122);imshow(II)

I = imread('E:\Underwater target recognition\chinamm2019uw_train\jpegimage\000000.jpg');
figure;imshow(I)

% extract RGB
R = I(:, :, 1);
G = I(:, :, 2);
B = I(:, :, 3);
R0 = double(R);
G0 = double(G);
B0 = double(B);

[N1, M1] = size(R);

Rlog = log(R0+1);
Rfft2 = fft2(R0);

sigma1 = 128;
F1 = fspecial('gaussian', [N1,M1], sigma1);
Efft1 = fft2(double(F1));

DR0 = Rfft2.* Efft1;
DR = ifft2(DR0);

DRlog = log(DR +1);
Rr1 = Rlog - DRlog;

sigma2 = 256;
F2 = fspecial('gaussian', [N1,M1], sigma2);
Efft2 = fft2(double(F2));

DR0 = Rfft2.* Efft2;
DR = ifft2(DR0);

DRlog = log(DR +1);
Rr2 = Rlog - DRlog;

sigma3 = 512;
F3 = fspecial('gaussian', [N1,M1], sigma3);
Efft3 = fft2(double(F3));

DR0 = Rfft2.* Efft3;
DR = ifft2(DR0);

DRlog = log(DR +1);
Rr3 = Rlog - DRlog;

Rr = (Rr1 + Rr2 +Rr3)/3;

a = 125;
II = imadd(R0, G0);
II = imadd(II, B0);
Ir = immultiply(R0, a);
C = imdivide(Ir, II);
C = log(C+1);

Rr = immultiply(C, Rr);
EXPRr = exp(Rr);
MIN = min(min(EXPRr));
MAX = max(max(EXPRr));
EXPRr = (EXPRr - MIN)/(MAX - MIN);
EXPRr = adapthisteq(EXPRr);

Glog = log(G0+1);
Gfft2 = fft2(G0);

DG0 = Gfft2.* Efft1;
DG = ifft2(DG0);

DGlog = log(DG +1);
Gg1 = Glog - DGlog;


DG0 = Gfft2.* Efft2;
DG = ifft2(DG0);

DGlog = log(DG +1);
Gg2 = Glog - DGlog;


DG0 = Gfft2.* Efft3;
DG = ifft2(DG0);

DGlog = log(DG +1);
Gg3 = Glog - DGlog;

Gg = (Gg1 + Gg2 +Gg3)/3;

Ig = immultiply(G0, a);
C = imdivide(Ig, II);
C = log(C+1);

Gg = immultiply(C, Gg);
EXPGg = exp(Gg);
MIN = min(min(EXPGg));
MAX = max(max(EXPGg));
EXPGg = (EXPGg - MIN)/(MAX - MIN);
EXPGg = adapthisteq(EXPGg);

Blog = log(B0+1); 
Bfft2 = fft2(B0); 

DB0 = Bfft2.* Efft1; 
DB = ifft2(DB0); 

DBlog = log(DB +1); 
Bb1 = Blog - DBlog; 

DB0 = Gfft2.* Efft2; 
DB = ifft2(DB0); 
DBlog = log(DB +1); 
Bb2 = Blog - DBlog;
DB0 = Gfft2.* Efft3; 
DB = ifft2(DB0);
DBlog = log(DB +1); 
Bb3 = Blog - DBlog;
Bb = (Bb1 + Bb2 + Bb3)/3; 

Ib = immultiply(B0, a);
C = imdivide(Ib, II); 
C = log(C+1);
Bb = immultiply(C, Bb); 
EXPBb= exp(Bb);
MIN = min(min(EXPBb)); 
MAX = max(max(EXPBb));
EXPBb = (EXPBb - MIN)/(MAX - MIN);
EXPBb = adapthisteq(EXPBb);

result = cat(3, EXPRr, EXPGg, EXPBb);  
figure;
subplot(121), imshow(I);
subplot(122), imshow(result);
figure;imshow(result)

Result:

Above are some column algorithm implementations of Retinex.

Published 15 original articles, won praise 5, visited 4341
Private letter follow

Tags: Attribute MATLAB

Posted on Tue, 11 Feb 2020 05:20:20 -0500 by curtisdw