Main Content

Deblur Images Using Regularized Filter

This example shows how to use regularized deconvolution to deblur images. Regularized deconvolution can be used effectively when limited information is known about the additive noise and constraints (such as smoothness) are applied on the recovered image. The blurred noisy image is restored by a constrained least square restoration algorithm that uses a regularized filter.

Simulate Gaussian Blur and Gaussian Noise

Read and display a pristine image that does not have blur or noise.

I = im2double(imread("tissue.png"));
imshow(I)
title("Original Image")
text(size(I,2),size(I,1)+15, ...
    "Image courtesy of Alan Partin, Johns Hopkins University", ...
    FontSize=7,HorizontalAlignment="right")

Figure contains an axes object. The axes object with title Original Image contains 2 objects of type image, text.

Simulate a blurred image that might result from an out-of-focus lens. First, create a point-spread function, PSF, by using the fspecial function and specifying a Gaussian filter of size 11-by-11 and standard deviation 5. Then, convolve the point-spread function with the image by using imfilter.

PSF = fspecial("gaussian",11,5);
blurred = imfilter(I,PSF,"conv");

Add zero-mean Gaussian noise to the blurred image by using the imnoise function.

noise_mean = 0;
noise_var = 0.02;
blurred_noisy = imnoise(blurred,"gaussian",noise_mean,noise_var);

Display the blurred noisy image.

imshow(blurred_noisy)
title("Blurred and Noisy Image")

Figure contains an axes object. The axes object with title Blurred and Noisy Image contains an object of type image.

Restore Image Using Estimated Noise Power

Restore the blurred image by using the deconvreg function, supplying the noise power (NP) as the third input parameter. To illustrate how sensitive the algorithm is to the value of noise power, this example performs three restorations.

For the first restoration, use the true NP. Note that the example outputs two parameters here. The first return value, reg1, is the restored image. The second return value, lagra, is a scalar Lagrange multiplier on which the regularized deconvolution has converged. This value is used later in the example.

NP = noise_var*numel(I);
[reg1,lagra] = deconvreg(blurred_noisy,PSF,NP);
imshow(reg1)
title("Restored with True NP")

Figure contains an axes object. The axes object with title Restored with True NP contains an object of type image.

For the second restoration, use a slightly overestimated noise power. The restoration has poor resolution.

NP_scale1 = 1.1;
reg2 = deconvreg(blurred_noisy,PSF,NP*NP_scale1);
imshow(reg2)
title("Restored with Larger NP")

Figure contains an axes object. The axes object with title Restored with Larger NP contains an object of type image.

For the third restoration, use a slightly underestimated noise power. The smaller the noise power multiplier, the stronger the noise amplification and ringing from the image borders.

NP_scale2 = 0.9;
reg3 = deconvreg(blurred_noisy,PSF,NP*NP_scale2);
imshow(reg3)
title("Restored with Smaller NP")

Figure contains an axes object. The axes object with title Restored with Smaller NP contains an object of type image.

Reduce Noise Amplification and Ringing

You can reduce the noise amplification and ringing along the boundary of the image by calling the edgetaper function prior to deconvolution. The image restoration becomes less sensitive to the noise power parameter.

tapered = edgetaper(blurred_noisy,PSF);
reg4 = deconvreg(tapered,PSF,NP*NP_scale2);
imshow(reg4)
title("Restored with Smaller NP and Edge Tapering")

Figure contains an axes object. The axes object with title Restored with Smaller NP and Edge Tapering contains an object of type image.

Use Lagrange Multiplier

Restore the blurred and noisy image, assuming that the optimal solution is already found and the corresponding Lagrange multiplier is known. In this case, any value passed for noise power, NP, is ignored.

To illustrate how sensitive the algorithm is to the Lagrange multiplier, this example performs three restorations. The first restoration uses the lagra output from the reg1 restoration performed earlier.

reg5 = deconvreg(tapered,PSF,[],lagra);
imshow(reg5)
title("Restored with LAGRA")

Figure contains an axes object. The axes object with title Restored with LAGRA contains an object of type image.

The second restoration uses a larger value than the lagra output from the reg1 restoration. A larger value increases the significance of the constraint. By default, this leads to oversmoothing of the image.

lagra_scale1 = 100;
reg6 = deconvreg(tapered,PSF,[],lagra*lagra_scale1);
imshow(reg6)
title("Restored with Large LAGRA")

Figure contains an axes object. The axes object with title Restored with Large LAGRA contains an object of type image.

The third restoration uses a smaller value than the lagra output from the reg1 restoration. A small value weakens the constraint (the smoothness requirement set for the image). It amplifies the noise. For the extreme case when the Lagrange multiplier equals 0, the reconstruction is a pure inverse filtering.

lagra_scale2 = 0.01;
reg7 = deconvreg(tapered,PSF,[],lagra*lagra_scale2);
imshow(reg7)
title("Restored with Small LAGRA")

Figure contains an axes object. The axes object with title Restored with Small LAGRA contains an object of type image.

Use Different Smoothness Constraint

Restore the blurred noisy image using a different constraint for the regularization operator. Instead of using the default Laplacian constraint on image smoothness, constrain the image smoothness using a 3-by-3 Laplacian of Gaussian (LoG) operator.

h = fspecial("log",3);
reg8 = deconvreg(blurred_noisy,PSF,[],lagra,h);
imshow(reg8)
title("Constrained by LoG Operator")

Figure contains an axes object. The axes object with title Constrained by LoG Operator contains an object of type image.

See Also

| | |

Related Topics