Lib4U

‎"Behind every stack of books there is a flood of knowledge."

Denoising Using Linear Filtering

denoise-2d-err-lena

Denoising of 1D signals.

    • Load a signal and add noise. We recommand that you test the 1D denoising on ‘Piece-Regular’ signal and on ‘gaussiannoise’. Download functionload_signal.

      n = 1024;
      name = ‘gaussiannoise’; % filtered gaussian noise
      name = ‘piece-regular’; % piecewise regular
      f = load_signal(name, n);
      sigma = 0.03 * (max(f)-min(f)); % noise level
      fn = f + sigma*randn(1,n); % noisy signal
      % plot signals
      subplot(2,1,1);
      plot(f); axis tight; title(‘Original’);
      subplot(2,1,2);
      plot(fn); axis tight; title(‘Noisy’);

    • Filtering with a gaussian kernel of variance mu (download functions build_gaussian_filter and perform_convolution)

      % mu is width of the filter, approximat
      mu = 5;
      m = 101; % total size of the filter, odd number is better
      h = build_gaussian_filter(m,mu/(4*n),n);
      fnh = perform_convolution(fn,h);
      % display
      subplot(3,1,1);
      plot(f); axis tight; title(‘Original’);
      subplot(3,1,2);
      plot(fn); axis tight; title(‘Noisy’);
      subplot(3,1,3);
      plot(fnh); axis tight; title(‘Estimated’);

       
      Denoising with linear filtering.You can see that blurring tends to destroy singularities on the left. On the contrary, linear filtering is well adapted to stationary signals such as filtered noise (right). Linear denoising is bound to fail on complex signals (see the next lecture for a non-linear denoising using wavelets).
    • Test various values for the variance and choose the best one.

mu_list = linspace( 0, 12, 15 );
err = [];
clf;
for i = 1:length(mu_list)
mu = mu_list(i);
h = build_gaussian_filter(m,mu/(4*n),n);
fnh = perform_convolution(fn,h);
% compute the error
e = sum( (f(:)-fnh(:)).^2 );
err = [ err, e ];
end
save_image([rep ‘denoise-progr-‘ name]);
[tmp,i] = min(err); mu = mu_list(i);
% plot the denoising result here

 
Dectection of the optimal value for the variance of the filter. This value is quite low, which result in a conservative filter. This filter prefers not to remove too much noise.
    • We can compare the result of the best gaussian filter with the Wiener filter, which achieve the best denoising using filtering. This asumes that the signal is drawn from a translation invariant (stationary) random vector. To denoise gaussian stationary signals, this wiener filter is also the best denoising one could achieved (including non-linear process). This is of course not true for piece-wise regular signals. Note that one should take special care at the boundary since FFT asume periodic condition. This is why we crop the image before compute the PSNR.

ff = fft(f); % Fourier transform of the signal
ffn = fft(fn); % Fourier transform of the noisy input
pf = abs(ff).^2; % spectral power
% fourier transform of the wiener filter
hwf = pf./(pf+ n*sigma^2);
% filter in the spatial domain (just for display)
hw = fftshift( ifft(hwf) );
% perform the filtering over the fourier domain
fw = real( ifft(ffn .* hwf) );
% compute the error. Warning: remove the boundary.
sel = n/2-256:n/2+256;
pwiener = psnr(f(sel),fw(sel));
pgauss = psnr(f(sel),fnh(sel));
% display

 
Comparison of Gaussian filtering vs. Wiener filtering.
Wiener is optimal among linear process.
But one could do better with non-linear denoising for
non gaussian sources.

Image denoising.

    • Load an image and add noise. Download files load_imagerescaleclamp and lena.

      n = 256;
      name = ‘gaussiannoise’;
      name = ‘lena’;
      I = load_image(name);
      % reduce size to speed up
      I = I(end/2-n/2+1:end/2+n/2,end/2-n/2+1:end/2+n/2);
      % to avoid saturation
      I = rescale(I,15,240);
      % add noise
      sigma = 0.18 * (max(I(:))-min(I(:)));
      In = I + sigma*randn(n);

    • Denoising with linear filtering.

      mu = 5; m = 31;
      h = build_gaussian_filter([m m],mu/(4*n),[n n]);
      Inh = perform_convolution(In,h); % calcule la convolution
      subplot(1,3,1);
      image( I ); axis image; axis off; title(‘Original’);
      subplot(1,3,2);
      image( clamp(In,0,255) ); axis image; axis off; title(‘Noisy’);
      subplot(1,3,3);
      image( clamp(Inh,0,255) ); axis image; axis off; title(‘Smoothed’);
      colormap gray(256);


      Denoising with linear filtering.
      Up: lena image.
      Bottom: filtered gaussian noise.
    • Test a large number of variance mu and keep the best.

mu_list = linspace( 0, 8, 12 );
err = [];
for mu = mu_list

end
% display the error function
plot(err, ‘.-‘)
% selection of the best variance
[tmp,i] = min(err); mu = mu_list(i);

display the denoising result

 
Progression of the denoising for increasing filtering size.

 

      
      
Denoising result for the optimal value of the variance of the filters.
Once again, there is still a lot of noise after the blurring.

 

    • One can also compare the best gaussian filter with the wiener filter.

fI = fft2(I);
fIn = fft2(In);
pf = abs(fI).^2; % spectral power
% fourier transform of the wiener filter
hwf = pf./(pf+ n^2*sigma^2);
% filter in the spatial domain
hw = fftshift( ifft2(hwf) );
% perform the filtering over the fourier
Iw = real( ifft2(fIn .* hwf) );
sel = n/2-100:n/2+100;
pwiener = psnr(I(sel,sel),Iw(sel,sel));
pgauss = psnr(I(sel,sel),Inh(sel,sel));
% display

 
Comparison of best gaussian filter with Wiener filter.
As in 1D, Wiener filter is the best linear filtering denoising.
Note however that non-linear denoisings achieve much better results on
natural images.

Source:

http://www.ceremade.dauphine.fr/~peyre/teaching/wavelets/tp3.html

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

Virtual Fashion Technology

Virtual Fashion Education

toitocuaanhem

"chúng tôi chỉ là tôi tớ của anh em, vì Đức Kitô" (2Cr 4,5b)

VentureBeat

News About Tech, Money and Innovation

digitalerr0r

Modern art using the GPU

Theme Showcase

Find the perfect theme for your blog.

lsuvietnam

Learn to Learn

Gocomay's Blog

Con tằm đến thác vẫn còn vương tơ

Toán cho Vật lý

Khoa Vật lý, Đại học Sư phạm Tp.HCM - ĐT :(08)-38352020 - 109

Maths 4 Physics & more...

Blog Toán Cao Cấp (M4Ps)

Bucket List Publications

Indulge- Travel, Adventure, & New Experiences

Lib4U

‎"Behind every stack of books there is a flood of knowledge."

The WordPress.com Blog

The latest news on WordPress.com and the WordPress community.

%d bloggers like this: