Wiener filter output doesn't look good enough-CodePudding

I have an image, I add some noise to it and try to denoise it using wiener filter:

Which in case of additive white noise and no blurring simplifies to:

Here is my code according to above formula, however, there is almost no difference with the input image.

import cv2
import numpy as np
img=cv2.imread('Images/P3.jpg',0)
freq2 = np.fft.fft2(img)
mean = 0
var = 100
sigma = var**0.5
gauss = np.random.normal(mean,sigma,np.shape(img))
courrupted=img gauss
freq2h = np.fft.fft2(gauss)
courrupted[courrupted<0]=0
courrupted[courrupted>255]=255
courrupted=courrupted.astype(np.uint8)

crfre=np.fft.fftshift(np.fft.fft2(courrupted))
sf=np.abs(crfre)**2

wiener=sf/(sf (100))
F_hat = crfre*wiener
f_hat = np.fft.ifft2( (F_hat))
restored = abs(f_hat)
normalizedImg=np.ones(img.shape)
normalizedImg = cv2.normalize(restored,  normalizedImg, 0, 255, cv2.NORM_MINMAX)

cv2.imwrite('output.jpg',normalizedImg)
cv2.imwrite('input.jpg',courrupted)

This is ground truth image:

This is input:

And this is output:

CodePudding user response：

The following works for me in Python/OpenCV/Numpy. As @Cris Luengo suggested, you need to test various values for the noise value, because the value you need may not be exactly your Gaussian variance input value.

Input:

import cv2
import numpy as np

# read image as grayscale
img = cv2.imread('pandas_noisy.jpg',0)

# take dft
dft = np.fft.fft2(img)

# get power spectral density of dft = square of magnitude
# where abs of complex number is the magnitude
pspec = (np.abs(dft))**2
print(np.amin(pspec))
print(np.amax(pspec))

# estimate noise power spectral density
# try different values to achieve compromise between noise reduction and softening/blurring
#noise = 100000000
#noise = 500000000
#noise = 1000000000
noise = 5000000000

# do wiener filtering
wiener = pspec/(pspec noise)
wiener = wiener*dft

# do dft to restore
restored = np.fft.ifft2(wiener)

# take real() component (or do abs())
restored = np.real(restored)
print(np.amin(restored))
print(np.amax(restored))

# clip and convert to uint8
restored = restored.clip(0,255).astype(np.uint8)

# save results
#cv2.imwrite('pandas_noisy_restored_100000000.jpg',restored)
#cv2.imwrite('pandas_noisy_restored_500000000.jpg',restored)
#cv2.imwrite('pandas_noisy_restored_1000000000.jpg',restored)
cv2.imwrite('pandas_noisy_restored_5000000000.jpg',restored)

# display results
cv2.imshow("input", img)
cv2.imshow("restored", restored)
cv2.waitKey(0)

Restored Result for Noise=100000000:

Restored Result for Noise=500000000:

Restored Result for Noise=1000000000:

Restored Result for Noise=5000000000: