3
$\begingroup$

I want to remove a noise for an image using MATLAB, when the observed image is $$f=u+v$$ where $u$ is the restored image (is the image i want recovered) and $v$ is the gaussian noise.

To restore $u$, I solve the following minimization problem: $$\min_{u \in H^1(\Omega)} \int_\Omega \gamma|\nabla u(x)|^2dx+ \int_\Omega (f(x)-u(x))^2dx,$$ where $\gamma$ is the regularization coefficient and $\Omega=[0,n]\times[0,m]$ and $[n,m]=$size($u$).

I want to solve the PDE (Euler-Lagrange) using MATLAB: \begin{eqnarray} div(\gamma \nabla u) + u = f \;in \;\Omega \\ \frac{\partial u}{\partial n}=0 \; in \; \partial \Omega \end{eqnarray}

Can anyone help me to solve this problem? Thank you!

I tried the following code :

clear all, close all,clc

uor=imread('gourd.bmp'); % the original image 
u0 = imnoise(uor,'gaussian',0,0.01);
u0=double(u0);
[m n]=size(u0);
uor=double(uor);
u=u0;
c=0.028;
h=1;     
for Iter=1:50, 
    for i=2:m-1,
      for j=2:n-1,
          Lap=0.003*(u(i+1,j)+u(i-1,j)-4*u(i,j)+u(i,j+1)+u(i,j-1));
          u(i,j)=(u0(i,j)+(1/(2*c*h*h))*Lap);
     end
    end
    for i=2:m-1,
          u(i,1)=u(i,2);
          u(i,n)=u(i,n-1);
        end

    for j=2:n-1,
          u(1,j)=u(2,j);
          u(m,j)=u(m-1,j);
        end

        u(1,1)=u(2,2);
        u(1,n)=u(2,n-1); 
        u(m,1)=u(m-1,2);
        u(m,n)=u(m-1,n-1);

en=0.0;  
    for i=2:m-1,
      for j=2:n-1,
      ux=(u(i+1,j)-u(i,j))/h;
      uy=(u(i,j+1)-u(i,j))/h;
      fidelity=(u0(i,j)-u(i,j))*(u0(i,j)-u(i,j));

      en=en+c*fidelity;
      end
    end


Energy(Iter)=en; 

%  Error between uor and u0
 ur=reshape(u,m*n,1);
 uori=reshape(uor,m*n,1);
 residu=norm(ur-uori)/norm(uori);

 [peaksnr, snr] = psnr(uor, u);

disp(['    iter ' num2str(Iter), ' :     Error = ' num2str(residu), ...
    ' ,    Peak-snr ' num2str(-peaksnr), ' ,    SNR ' num2str(snr)]);


  end 

% show the structural similarity index for measuring image quality 
[ssimval, mapssim] = ssim(u,uor);
disp([' the structural similarity index is ' num2str(ssimval)]);
figure,imshow(mapssim,[]); axis square; 

figure,imagesc(u); axis image; axis off; colormap(gray);

The original image is here : https://www.dropbox.com/s/4bccby1f4lxp4j9/gourd.rar?dl=0

gourd

Best regards

$\endgroup$
1
$\begingroup$

You're trying to solve what's called Perona Malik Non Linear Diffusion Problem (Sometimes people call it, by mistake, Anisotropic Diffusion).

Anyhow, the simplest code for that is Anisotropic Diffusion (Perona & Malik) on The MATLAB File Exchange.

There is a more advanced (Anisotropic for real) algorithm in Fast Anisotropic Curvature Preserving Smoothing (Also on the File Exchange - Fast Anisotropic Curvature Preserving Smoothing).

Enjoy...

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.