2

我知道，卷积original_image * filter = blur_image在哪里。*因此，filter = ifft(fft(blur)/fft(original))

我有一个原始图像、已知的过滤器和已知的模糊图像。我尝试了以下代码。我只想比较使用 fft 和 ifft 计算的过滤器并将其与已知过滤器进行比较。我在 Matlab 中试过：

orig = imread("orig.png")
blur = imread("blur.png")
fftorig = fft(orig)
fftblur = fft(blur)
div = fftblur/fftorig
conv = ifft(div)

结果没有意义。我看到它div包含许多NaN值，fftblur并且fftorig都包含许多 0 值。我需要对此做些什么吗？比如使用fftshift？

编辑：为了更容易理解，我现在使用的图片来自： http: //matlabgeeks.com/tips-tutorials/how-to-blur-an-image-with-a-fourier-transform-in- matlab-第i部分/

我决定使用以下链接计算origimageand的内核blurimageunpad：

kernelc = real(ifft2(fft2(origimage)./fft2(blurimageunpad));
imagesc(kernelc)
colormap gray

结果如下：

https://imgur.com/a/b7uvj

这显然与该链接顶部提到的高斯模糊不匹配

4

1 回答 1

4

这就是维纳滤波器派上用场的地方。它通常用于反卷积——从过滤后的图像和卷​​积核中估计原始的、未过滤的图像。但是，由于卷积的交换性，反卷积与您要解决的问题完全相同。也就是说，如果g = f * h，那么您可以从g和h（反卷积）估计f ，或者您可以从g和f估计h，以完全相同的方式。

反卷积

Wiener filter Wikipedia 页面在数学上很重要，但是以一种简单的方式，您会寻找内核具有较小值的频率（与图像中的噪声相比），并且您不会在这些频率上应用除法。这称为正则化，您可以在其中施加一些约束以避免噪声主导结果。

这是 MATLAB 代码，in用于模糊输入图像和psf空间域滤波器：

% Create a PSF and a blurred image:
original = imread('cameraman.tif');
sz = size(original);
psf = zeros(sz);
psf(fix(sz(1)/2)+(-5:5),fix(sz(1)/2)+(-10:10)) = 1;
psf = psf/sum(psf(:));
% in = conv2(original,psf,'same'); % gives boundary issues, less of a problem for smaller psf
in = uint8(ifft2(fft2(original).*fft2(ifftshift(psf))));

% Input parameter:
K = 1e-2; 

% Compute frequency-domain PSF:
H = fft2(ifftshift(psf));

% Compute the Wiener filter:
cH = conj(H);
HcH = H .* cH;
K = K * max(max(abs(HcH)));
w = cH ./ (HcH + K);

% Apply the Wiener filter in the Frequency domain:
out = real(ifft2(w .* fft2(in)));

in这是三个不同值的维纳滤波器的图像和输出K：

可以看到，选对K是很重要的。如果它太低，则没有足够的正则化。如果它太高，那么就会有太多的伪影（去卷积不足）。这个参数K取决于输入图像，尽管有一些技巧可以估计它。此外，像这样的简单过滤器永远无法完美地消除我在这里放置的刺眼模糊。为了获得更好的反卷积，需要更先进的迭代方法。

估计内核

让我们反过来，psf从original和估计in。可以直接进行除法（相当于带有 的维纳滤波器K=0），但输出噪声很大。在原始图像的频域值非常低的情况下，您会得到很差的估计。选择正确K可以更好地估计内核。如果 aK太大，结果又是一个很差的近似值。

% Direct estimation by division
psf1 = fftshift(ifft2(fft2(in)./fft2(original)));

% Wiener filter approach
K = 1e-7;
H = fft2(original);
cH = conj(H);
HcH = H .* cH;
K = K * max(max(abs(HcH)));
w = cH ./ (HcH + K);
psf2 = fftshift(real(ifft2(w .* fft2(in))));

（放大观察噪音）

编辑

我从您链接的网站下载了图像。我使用第一个结果，没有填充，并尽可能地从图像中裁剪帧，只留下数据和准确的数据：

original = imread('original.png');
original = original(33:374,120:460,1);

in = imread('blur_nopad.png');
in = in(33:374,120:460,1);

然后我完全按照上面的代码使用了代码，K一切都一样，得到了一个很好的结果，显示了一个稍微偏移的高斯核。

然后我对第二个结果（填充后）重复了同样的操作，得到了一个更差的结果，但仍然非常容易识别为高斯。

于 2018-02-16T07:28:09.740 回答