This example shows how to use Wiener deconvolution to deblur images. Wiener deconvolution can be useful when the point-spread function and noise level are either known or estimated.
|On this page…|
I = im2double(imread('cameraman.tif')); imshow(I); title('Original Image (courtesy of MIT)');
Simulate a blurred image that you might get from camera motion. Create a point-spread function, PSF, corresponding to the linear motion across 31 pixels (LEN=31), at an angle of 11 degrees (THETA=11). To simulate the blur, convolve the filter with the image using imfilter.
LEN = 21; THETA = 11; PSF = fspecial('motion', LEN, THETA); blurred = imfilter(I, PSF, 'conv', 'circular'); imshow(blurred); title('Blurred Image');
The simplest syntax for deconvwnr is deconvwnr(A, PSF, NSR), where A is the blurred image, PSF is the point-spread function, and NSR is the noise-power-to-signal-power ratio. The blurred image formed in Step 2 has no noise, so we'll use 0 for NSR.
wnr1 = deconvwnr(blurred, PSF, 0); imshow(wnr1); title('Restored Image');
Now let's try adding noise.
noise_mean = 0; noise_var = 0.0001; blurred_noisy = imnoise(blurred, 'gaussian', ... noise_mean, noise_var); imshow(blurred_noisy) title('Simulate Blur and Noise')
In our first restoration attempt, we'll tell deconvwnr that there is no noise (NSR = 0). When NSR = 0, the Wiener restoration filter is equivalent to an ideal inverse filter. The ideal inverse filter can be extremely sensitive to noise in the input image, as the next image shows:
wnr2 = deconvwnr(blurred_noisy, PSF, 0); imshow(wnr2) title('Restoration of Blurred, Noisy Image Using NSR = 0')
The noise was amplified by the inverse filter to such a degree that only the barest hint of the man's shape is visible.
In our second attempt we supply an estimate of the noise-power-to-signal-power ratio.
signal_var = var(I(:)); wnr3 = deconvwnr(blurred_noisy, PSF, noise_var / signal_var); imshow(wnr3) title('Restoration of Blurred, Noisy Image Using Estimated NSR');
Even a visually imperceptible amount of noise can affect the result. Let's try keeping the input image in uint8 representation instead of converting it to double.
I = imread('cameraman.tif'); class(I)
ans = uint8
If you pass a uint8 image to imfilter, it will quantize the output in order to return another uint8 image.
blurred_quantized = imfilter(I, PSF, 'conv', 'circular'); class(blurred_quantized)
ans = uint8
Again, we'll try first telling deconvwnr that there is no noise.
wnr4 = deconvwnr(blurred_quantized, PSF, 0); imshow(wnr4) title('Restoration of blurred, quantized image using NSR = 0');
Next, we supply an NSR estimate to deconvwnr.
uniform_quantization_var = (1/256)^2 / 12; signal_var = var(im2double(I(:))); wnr5 = deconvwnr(blurred_quantized, PSF, ... uniform_quantization_var / signal_var); imshow(wnr5) title('Restoration of Blurred, Quantized Image Using Computed NSR');