9

I am trying to implement the Wiener Filter to perform deconvolution on blurred image. My implementation is like this

import numpy as np
from numpy.fft import fft2, ifft2

def wiener_filter(img, kernel, K = 10):
    dummy = np.copy(img)
    kernel = np.pad(kernel, [(0, dummy.shape[0] - kernel.shape[0]), (0, dummy.shape[1] - kernel.shape[1])], 'constant')
    # Fourier Transform
    dummy = fft2(dummy)
    kernel = fft2(kernel)
    kernel = np.conj(kernel) / (np.abs(kernel) ** 2 + K)
    dummy = dummy * kernel
    dummy = np.abs(ifft2(dummy))
    return np.uint8(dummy)

This implementation is based on the Wiki Page.

The TIFF image used is from : http://www.ece.rice.edu/~wakin/images/lena512color.tiff
But here is a PNG version:

I have a input image motion blurred by a diagonal kernel and some gaussian additive noise is added to it. The lena picture is 512x512 and the blurring kernel is 11x11.

When I apply my wiener_filter to this image the result is like this. enter image description here.

I think this deblurred image is not of good quality. So I would like to ask if my implementation is correct.

Update the way I add noise.

from scipy.signal import gaussian, convolve2d

def blur(img, mode = 'box', block_size = 3):
    # mode = 'box' or 'gaussian' or 'motion'
    dummy = np.copy(img)
    if mode == 'box':
        h = np.ones((block_size, block_size)) / block_size ** 2
    elif mode == 'gaussian':
        h = gaussian(block_size, block_size / 3).reshape(block_size, 1)
        h = np.dot(h, h.transpose())
        h /= np.sum(h)
    elif mode == 'motion':
        h = np.eye(block_size) / block_size
    dummy = convolve2d(dummy, h, mode = 'valid')
    return np.uint8(dummy), h

def gaussian_add(img, sigma = 5):
    dummy = np.copy(img).astype(float)
    gauss = np.random.normal(0, sigma, np.shape(img))
    # Additive Noise
    dummy = np.round(gauss + dummy)
    # Saturate lower bound
    dummy[np.where(dummy < 0)] = 0
    # Saturate upper bound
    dummy[np.where(dummy > 255)] = 255
    return np.uint8(dummy)
Cris Luengo
  • 55,762
  • 10
  • 62
  • 120
yc2986
  • 1,101
  • 3
  • 20
  • 23
  • just out of interest.. what happens when you run a "sharpen" convolution on the new (blurred) image. Does it fix? – VC.One Feb 04 '16 at 05:01
  • I am new to image processing. Is the sharpen kernel you mean here something like (lapacian matrix + 1 centered by 0s matrix)? Thanks. @VC.One – yc2986 Feb 04 '16 at 06:07
  • Nevermind I tested and it just sharpens the (angled) motion blur so it's more crisp but also looks "dirty". – VC.One Feb 04 '16 at 06:37

3 Answers3

3

Use skimage.restoration.wiener, which is usually used like:

>>> from skimage import color, data, restoration
>>> img = color.rgb2gray(data.astronaut())
>>> from scipy.signal import convolve2d
>>> psf = np.ones((5, 5)) / 25
>>> img = convolve2d(img, psf, 'same')
>>> img += 0.1 * img.std() * np.random.standard_normal(img.shape)
>>> deconvolved_img = restoration.wiener(img, psf, 1100)

I have also used it in: Deblur an image using scikit-image.

Community
  • 1
  • 1
gsamaras
  • 71,951
  • 46
  • 188
  • 305
2

For data comparison, you can find a sample implementation of Wiener filtering and unsupervisived Wiener filtering at

http://scikit-image.org/docs/dev/auto_examples/plot_restoration.html

If you give your original image data, we may be able to help further.

EDIT: Original link seems to be down, try this one: http://scikit-image.org/docs/dev/auto_examples/filters/plot_restoration.html

tfv
  • 6,016
  • 4
  • 36
  • 67
  • My input is a colorful lena image provided by this link: http://www.ece.rice.edu/~wakin/images/lena512color.tiff. I was loading this image using PIL Image and convert it to RGB. I added a convolutional noise and another additive gaussian noise with a relatively small standard deviation. – yc2986 Feb 04 '16 at 05:51
  • 2
    @tfv: Your link is down – CpCd0y Mar 20 '17 at 08:48
2

We could try unsupervised weiner too (deconvolution with a Wiener-Hunt approach, where the hyperparameters are automatically estimated, using a stochastic iterative process (Gibbs sampler), as described here):

deconvolved, _ = restoration.unsupervised_wiener(im, psf)
Sandipan Dey
  • 21,482
  • 2
  • 51
  • 63