[Diffusion filtering for image denoising] Image linear diffusion filtering, edge enhancement linear and nonlinear anisotropic filtering (Matlab code implementation)

Welcome to this blog

Advantages of bloggers:Blog content should be as thoughtful and logical as possible for the convenience of readers.

Motto:He who travels a hundred miles is half as good as ninety.

The directory of this article is as follows:

Table of Contents

1 Overview

2 Operation results

2.1 non_aniso

2.2 inhomo_iso

2.3 heat_imp

2.4 heat_explicit

3 References

4 Matlab code implementation


1 Overview

[Diffusion filtering for image denoising] Image linear diffusion filtering, edge enhancement linear and nonlinear anisotropic filtering

Various diffusion-based image filtering methods:
1. Linear diffusion filtering using thermal equations – solved using implicit and explicit Euler methods.
2. Edge-enhanced linear anisotropic diffusion filtering.
3. Edge-enhanced nonlinear anisotropic diffusion filtering.

In the field of image processing, removing noise and improving details in images is a very important task. To this end, various diffusion-based image filtering methods have been developed, aiming to improve image quality and enhance image features. Here are some common filtering methods:

1. Linear diffusion filtering: Use the thermal equation model to perform linear diffusion filtering. This method can be solved using implicit and explicit Euler methods to improve image smoothness and reduce noise. Linear diffusion filtering can achieve image denoising and smoothing by adaptively adjusting the brightness values of neighboring pixels.

2. Edge-enhanced linear anisotropic diffusion filtering: This method combines the ideas of linear diffusion filtering and anisotropic filtering to highlight the edge features in the image while maintaining image smoothness. By weighting the gradients of neighborhood pixels, edge-enhanced linear anisotropic diffusion filtering can effectively improve the edge details of the image, making the image clearer and sharper.

3. Edge-enhanced nonlinear anisotropic diffusion filtering: Compared with linear filtering, nonlinear filtering methods are more flexible and accurate when processing images. Edge-enhanced nonlinear anisotropic diffusion filtering makes the processed image more advantageous in terms of detail and edge preservation by considering the nonlinear relationship between pixels. This method can remove noise while maintaining the details and texture information of the image.

These diffusion-based image filtering methods provide powerful tools for image denoising and enhancement. According to different application requirements and image characteristics, choosing an appropriate filtering method can significantly improve the image quality and enhance the visualization effect of the image. Whether in the fields of scientific research, medical imaging or graphic design, these filtering methods play an important role in bringing users a better image processing experience.

2 Operation results

2.1 non_aniso

2.2 inhomo_iso

2.3 heat_imp

2.4 heat_explicit

Part of the code:

% set up finite difference parameters
alpha =.5;
k = 1;
h = 1;

lambda = (alpha^2)*(k/(h^2));

[m n] = size(w);
%stick image into a vector
w_vec = reshape(w,n*n,1);
w_old = w;
w_new = w;

%smooth the image
im_smth = filter_function(w,1);
im_smth = im_smth';
% required for g() calculation
[dx_im_smth dy_im_smth] = gradient(im_smth);
gr_im_smth = dx_im_smth.^2 + dy_im_smth.^2;

[mmm nnn] = size(gr_im_smth);
% g() calculation
for i=1:mmm
    for j=1:nnn
        g(i,j) = 1/(1 + (gr_im_smth(i,j)/32));
    end
end
%g = g';
jj=1;

[im_sx im_sy] = gradient(im_smth);
% diffusion tensor D is precomputed here
D = zeros(n*n,2,2);
for i=1:n*n
    row = ceil(i/n);
    col = i - (row-1) * n;
  % if((col > 1) & amp; & amp; (col < n) & amp; & amp; (row > 1) & amp; & amp; (row < n))
        
        % choose eigenvector parallel and perpendicular to gradient
        eigen_vec = [im_sx(row,col) im_sy(row,col); im_sy(row,col) -im_sx(row,col) ];
        %chooseeigenvalues
        eigen_vec(:,1) = eigen_vec(:,1) ./ norm(eigen_vec(:,1));
        eigen_vec(:,2) = eigen_vec(:,2) ./ norm(eigen_vec(:,2));
        eigen_val = [g(row,col) 0;0 1];
        
        % form diffusion tensor
        D(i,:,:) = eigen_vec * eigen_val * (eigen_vec');
        % end
end

figure;
for k=1:400 % for each iteration
    for i=1:n*n % solve using Jacobi iterations scheme
        row = ceil(i/n); %compute what row this pixel belongs to in original image
        col = i - (row-1) * n; % compute cols similarly

        %different if conditions handles pixels at different locations in
        %the image as depending on their location they may or may not have
        %all their neighbor pixels which will be required for finite
        %differences
        
        if((col > 1) & amp; & amp; (col < n) & amp; & amp; (row > 1) & amp; & amp; (row < n))
            
            s = -lambda * ((D(i,1,1) * w_old(i-1)) + ((D(i,1,1) - D(i,1,2) - D(i,2, 1)) * (w_old(i + 1))) + ...
                (D(i,2,2) * w_old(i-n)) + ((D(i,2,2) - D(i,1,2) - D(i,2,1)) * (w_old(i + n))) + ((D(i,1,2) + D(i,2,1)) * w_old(i + n + 1)));
            w_new(i) = (-s + w_old(i))/(1 + lambda * (2 * D(i,1,1) + 2 * D(i,2,2) - D(i,1,2 ) - D(i,2,1)));
            
        elseif((row == 1) & amp; & amp; (col > 1) & amp; & amp; (col < n))
            s = -lambda * ((D(i,1,1) * w_old(i-1)) + ((D(i,1,1) - D(i,1,2) - D(i,2, 1)) * (w_old(i + 1))) + ...
                ((D(i,2,2) - D(i,1,2) - D(i,2,1)) * (w_old(i + n))) + ((D(i,1,2) + D(i,2,1)) * w_old(i + n + 1)));
            
            %s = -(lambda) * (w_old(i + 1) + w_old(i-1) + w_old(i + n));
            w_new(i) = (-s + w_old(i))/(1 + 4*lambda);
        elseif((row == n) & amp; & amp; (col > 1) & amp; & amp; (col < n))
            s = -lambda * ((D(i,1,1) * w_old(i-1)) + ((D(i,1,1) - D(i,1,2) - D(i,2, 1)) * (w_old(i + 1))) + ...
                (D(i,2,2) * w_old(i-n)));
            
            %s = -(lambda) * (w_old(i + 1) + w_old(i-1) + w_old(i-n));
            w_new(i) = (-s + w_old(i))/(1 + 4*lambda);
        elseif((col == 1) & amp; & amp; (row > 1) & amp; & amp; (row < n))
            s = -lambda * (((D(i,1,1) - D(i,1,2) - D(i,2,1)) * (w_old(i + 1))) + ...
                (D(i,2,2) * w_old(i-n)) + ((D(i,2,2) - D(i,1,2) - D(i,2,1)) * (w_old(i + n))) + ((D(i,1,2) + D(i,2,1)) * w_old(i + n + 1)));
            
            %s = -(lambda) * (w_old(i + 1) + w_old(i + n) + w_old(i-n));
            w_new(i) = (-s + w_old(i))/(1 + 4*lambda);
        elseif((col == n) & amp; & amp; (row > 1) & amp; & amp; (row < n))
            s = -lambda * ((D(i,1,1) * w_old(i-1)) + ...
                (D(i,2,2) * w_old(i-n)) + ((D(i,2,2) - D(i,1,2) - D(i,2,1)) * (w_old(i + n))));
            %s = -(lambda) * (w_old(i-1) + w_old(i + n) + w_old(i-n));
            w_new(i) = (-s + w_old(i))/(1 + 4*lambda);
        elseif((col==1) & amp; & amp; (row==1))

% set up finite difference parameters
alpha =.5;
k = 1;
h = 1;

lambda = (alpha^2)*(k/(h^2));

[m n] = size(w);
%stick image into a vector
w_vec = reshape(w,n*n,1);
w_old = w;
w_new = w;

%smooth the image
im_smth = filter_function(w,1);
im_smth = im_smth’;
% required for g() calculation
[dx_im_smth dy_im_smth] = gradient(im_smth);
gr_im_smth = dx_im_smth.^2 + dy_im_smth.^2;

[mmm nnn] = size(gr_im_smth);
% g() calculation
for i=1:mmm
for j=1:nnn
g(i,j) = 1/(1 + (gr_im_smth(i,j)/32));
end
end
%g = g’;
jj=1;

[im_sx im_sy] = gradient(im_smth);
% diffusion tensor D is precomputed here
D = zeros(n*n,2,2);
for i=1:n*n
row = ceil(i/n);
col = i – (row-1) * n;
% if((col > 1) & amp; & amp; (col < n) & amp; & amp; (row > 1) & amp; & amp; (row < n)) % choose eigenvector parallel and perpendicular to gradient
eigen_vec = [im_sx(row,col) im_sy(row,col); im_sy(row,col) -im_sx(row,col) ];
%chooseeigenvalues
eigen_vec(:,1) = eigen_vec(:,1) ./ norm(eigen_vec(:,1));
eigen_vec(:,2) = eigen_vec(:,2) ./ norm(eigen_vec(:,2));
eigen_val = [g(row,col) 0;0 1];

% form diffusion tensor
D(i,:,:) = eigen_vec * eigen_val * (eigen_vec’);
% end
end

figure;
for k=1:400 % for each iteration
for i=1:n*n % solve using Jacobi iterations scheme
row = ceil(i/n); %compute what row this pixel belongs to in original image
col = i – (row-1) * n; % compute cols similarly

%different if conditions handles pixels at different locations in
%the image as depending on their location they may or may not have
%all their neighbor pixels which will be required for finite
%differences

if((col > 1) & amp; & amp; (col < n) & amp; & amp; (row > 1) & amp; & amp; (row < n)) s = -lambda * ((D(i,1,1) * w_old(i-1)) + ((D(i,1,1) – D(i,1,2) – D(i,2, 1)) * (w_old(i + 1))) + …
(D(i,2,2) * w_old(i-n)) + ((D(i,2,2) – D(i,1,2) – D(i,2,1)) * (w_old(i + n))) + ((D(i,1,2) + D(i,2,1)) * w_old(i + n + 1)));
w_new(i) = (-s + w_old(i))/(1 + lambda * (2 * D(i,1,1) + 2 * D(i,2,2) – D(i,1,2 ) – D(i,2,1)));

elseif((row == 1) & amp; & amp; (col > 1) & amp; & amp; (col < n))
s = -lambda * ((D(i,1,1) * w_old(i-1)) + ((D(i,1,1) – D(i,1,2) – D(i,2, 1)) * (w_old(i + 1))) + …
((D(i,2,2) – D(i,1,2) – D(i,2,1)) * (w_old(i + n))) + ((D(i,1,2) + D(i,2,1)) * w_old(i + n + 1)));

%s = -(lambda) * (w_old(i + 1) + w_old(i-1) + w_old(i + n));
w_new(i) = (-s + w_old(i))/(1 + 4*lambda);
elseif((row == n) & amp; & amp; (col > 1) & amp; & amp; (col < n))
s = -lambda * ((D(i,1,1) * w_old(i-1)) + ((D(i,1,1) – D(i,1,2) – D(i,2, 1)) * (w_old(i + 1))) + …
(D(i,2,2) * w_old(i-n)));

%s = -(lambda) * (w_old(i + 1) + w_old(i-1) + w_old(i-n));
w_new(i) = (-s + w_old(i))/(1 + 4*lambda);
elseif((col == 1) & amp; & amp; (row > 1) & amp; & amp; (row < n))
s = -lambda * (((D(i,1,1) – D(i,1,2) – D(i,2,1)) * (w_old(i + 1))) + …
(D(i,2,2) * w_old(i-n)) + ((D(i,2,2) – D(i,1,2) – D(i,2,1)) * (w_old(i + n))) + ((D(i,1,2) + D(i,2,1)) * w_old(i + n + 1)));

%s = -(lambda) * (w_old(i + 1) + w_old(i + n) + w_old(i-n));
w_new(i) = (-s + w_old(i))/(1 + 4*lambda);
elseif((col == n) & amp; & amp; (row > 1) & amp; & amp; (row < n))
s = -lambda * ((D(i,1,1) * w_old(i-1)) + …
(D(i,2,2) * w_old(i-n)) + ((D(i,2,2) – D(i,1,2) – D(i,2,1)) * (w_old(i + n))));
%s = -(lambda) * (w_old(i-1) + w_old(i + n) + w_old(i-n));
w_new(i) = (-s + w_old(i))/(1 + 4*lambda);
elseif((col==1) & amp; & amp; (row==1))

3 References

Some content in the article is quoted from the Internet, and the source will be indicated or cited as a reference. It is inevitable that there will be some unfinished information. If there is anything inappropriate, please feel free to contact us to delete it.

[1] Zhang Jianwei, Wang Yihe, Chen Yunjie. Image denoising method based on nonlinear diffusion filtering structural information [J]. Computer Engineering and Design, 2016, 37(11):8.DOI:10.16208/j.issn1000-7024.2016 .11.021.

[2] Wen Wu, Miao Fang. Application of complex domain nonlinear diffusion filtering in image processing [J]. Microelectronics and Computers, 2012. DOI: CNKI: SUN: WXYJ.0.2012-06-015.

[3] Fu Yanli, Li Chao, Chen Hao, et al. Anisotropic diffusion filtering long-range acoustic logging image denoising method [J]. Applied Acoustics, 2022, 41(4):10.

[4]Xu Tao. Research on nonlinear diffusion image hybrid filtering denoising method[J]. Computer Simulation, 2020, 037(012):460-464.

[5] Li Zhiwei, Feng Xiangchu. Image denoising combining Wiener filtering and nonlinear diffusion [J]. Electronic Technology, 2007(9):4.DOI:10.3969/j.issn.1007-7820.2007.09.016.

4 Matlab code implementation