在 Matlab 问题中使用 FFT 进行二维反卷积答案

【问题标题】：2D Deconvolution using FFT in Matlab Problems在 Matlab 问题中使用 FFT 进行二维反卷积
【发布时间】：2013-10-09 15:20:39
【问题描述】：

我已经将我在 matlab 中创建的图像与我也在 matlab 中定义的 2D 高斯函数进行了卷积，现在我正在尝试对结果矩阵进行反卷积，以查看是否使用 fft2 和 ifft2 命令获得了 2D 高斯函数.然而，我得到的矩阵是不正确的（据我所知）。这是我到目前为止所做的代码：

% 输入图像代码 (img) [300x300 数组]

N = 100;
t = linspace(0,2*pi,50);
r = (N-10)/2;
circle = poly2mask(r*cos(t)+N/2+0.5, r*sin(t)+N/2+0.5,N,N);
img = repmat(circle,3,3);

% c = 0 sig = 1/64 (Z) [300x300 数组] 的二维高斯函数代码

x = linspace(-3,3,300);
y = x';
[X Y] = meshgrid(x,y);
Z = exp(-((X.^2)+(Y.^2))/(2*1/64));

% img 与 Z (C) [599x599 数组] 的二维卷积代码

C = conv2(img,Z);

% 我已经使用 img 和 C 的横截面轮廓向量测试了这个卷积是正确的，得到的 x-y 图是我对卷积的期望。

% 根据我对卷积的了解，该算法在傅立叶空间中用作乘数，因此通过将输出（卷积图像）的傅立叶变换除以输入（img），我应该得到点扩散函数（Z - 2D 高斯函数）经过傅里叶逆变换后通过除法应用于此结果。

% 用于尝试 2D 反卷积的代码

Fimg = fft2(img,599,599);

% 零填充添加以将结果增加到 599x599 数组

FC = fft2(C);
R = FC/Fimg;

% 我现在得到这个错误提示：警告：矩阵接近奇异或缩放错误。结果可能不准确。 RCOND = 2.551432e-22

iFR = ifft2(R);

我预计 iFR 会接近 Z，但我得到的结果完全不同。它可能是具有复数值的 Z 的近似值，但我似乎无法检查它，因为我不知道如何在 matlab 中绘制 3D 复矩阵。那么，如果有人能告诉我我的答案是正确还是不正确，以及如何让这种反卷积起作用？我将不胜感激。

【问题讨论】：

有点跑题了：你对信号进行卷积，而不是卷积
既然除法得到一个599x599的矩阵，那么如何得到300x300的二维高斯函数呢？
@IssamLaradji 结果图像为 599x599，2D 高斯包含在前 300x300 行中。我使用矩阵索引来取回一个 300x300 数组

标签： arrays matlab fft convolution ifft

【解决方案1】：

R = FC/Fimg 必须是 R = FC./Fimg; 您需要按元素进行除法。

【讨论】：

当然，这更有意义；多么菜鸟的错误。有没有办法将得到的复数矩阵绘制为曲面图，以便我可以将其与我的高斯进行比较？
@SeanJamesJamieson 我认为不应该有任何复杂的组件......不过不确定。但是以防万一您可以在绘制数组之前用 abs 包装数组。
是的，我这样做了，我为前 300 行和前 300 列的矩阵建立了索引，并以 4 x 10^-13 的最大逐点定义不确定性恢复了我的函数。谢谢大家的帮助
很好，但你怎么知道前 300 行和前 300 列代表高斯函数？例如，如果你取最后 300 行和列，会发生什么？我正在做类似的事情，这部分让我感到困惑：/
@IssamLaradji 恢复的高斯包含在反卷积后结果数组的第一象限内。我认为这与 matlab 填充初始数组的默认方式有关，因为现在当我在 padarray 函数中使用向量输入在数组的行和列之前和之后添加相等特定数量的零时，高斯包含在数组的中心

【解决方案2】：

这是该反卷积高斯的一些 Octave（版本 3.6.2）图。

% deconvolve in frequency domain
Fimg = fft2(img,599,599);
FC = fft2(C);
R = FC ./ Fimg;
r = ifft2(R);

% show deconvolved Gaussian
figure(1);
subplot(2,3,1), imshow(img), title('image');
subplot(2,3,2), imshow(Z), title('Gaussian');
subplot(2,3,3), imshow(C), title('image blurred by Gaussian');
subplot(2,3,4), mesh(X,Y,Z), title('initial Gaussian');
subplot(2,3,5), imagesc(real(r(1:300,1:300))), colormap gray, title('deconvolved Gaussian');
subplot(2,3,6), mesh(X,Y,real(r(1:300,1:300))), title('deconvolved Gaussian');

% show difference between Gaussian and deconvolved Gaussian
figure(2);
gdiff = Z - real(r(1:300,1:300));
imagesc(gdiff), colorbar, colormap gray, title('difference between initial Gaussian and deconvolved Guassian');

【讨论】：