MATLAB如何在频域实现Ram-Lak滤波器(Ramp filter)？

Question

我有一个任务要实现 Ram-Lak 滤波器，但几乎没有给出任何相关信息(除了查看 fft、ifft、fftshift、ifftshift)。

我有一个正弦图，我必须通过 Ram-Lak 过滤。还给出了投影的数量。

我尝试使用过滤器

                      1/4              if I == 0

(b^2)/(2*pi^2)  *     0                if I even

                      -1/(pi^2 * I^2)  if I odd

b好像是截止频率，我跟采样率有关系？

也有人说两个函数的卷积是傅立叶空间中的简单乘法。

我完全不明白如何实现过滤器，尤其是没有给出 b，没有告诉我我是什么，也不知道如何将它应用到正弦图中，我希望有人能在这里帮助我。我花了 2 个小时在谷歌上搜索并试图了解这里需要做什么，但我不明白如何实现它。

Answer 1

如果您想在傅立叶域中不进行过滤而进行逆拉东变换，则您列出的公式是一个中间结果。另一种方法是在空间域中使用卷积来完成整个滤波反投影算法，这在 40 年前可能更快；您最终会重新获得您发布的公式。但是，我现在不推荐它，尤其是对于您的第一次重建；您首先应该真正了解希尔伯特变换。

无论如何，这是一些执行强制性 Shepp-Logan 幻影滤波反投影重建的 Matlab 代码。我展示了如何使用 Ram-Lak 过滤器进行自己的过滤。如果我真的有动力，我会用一些 interp2 命令和求和来替换 radon/iradon。


phantomData=phantom();

N=size(phantomData,1);

θ = 0:179；
N_theta = 长度(θ)；
[R,xp] = radon(phantomData,theta);

% 做一个 Ram-Lak 过滤器。它只是 abs(f)。
N1 = 长度(xp);
freqs=linspace(-1, 1, N1).';
myFilter = abs(频率);
myFilter = repmat(myFilter, [1 N_theta]);

% 做我自己的 FT 域过滤
ft_R = fftshift(fft(R,[],1),1);
filteredProj = ft_R .* myFilter;
filteredProj = ifftshift(filteredProj,1);
ift_R = real(ifft(filteredProj,[],1));

% 告诉 matlab 在没有过滤器的情况下进行反 FBP
I1 = iradon(ift_R, theta, '线性', '无', 1.0, N);

子图(1,3,1)；imagesc(真实(I1))； title('手动过滤')
颜色图(灰色(256))；轴图像；轴关闭

% 为了比较，要求 matlab 使用他们的 Ram-Lak 过滤器实现
I2 = iradon(R, theta, 'linear', 'Ram-Lak', 1.0, N);

子图(1,3,2)；imagesc(真实(I2))； title('Matlab 滤波')
颜色图(灰色(256))；轴图像；轴关闭

% 为了好玩，故意重做过滤错误
% 排除高频以创建低分辨率重建
myFilter( myFilter > 0.1 ) = 0;
ift_R = real(ifft(ifftshift(ft_R .* myFilter,1),[],1));
I3 = iradon(ift_R, theta, '线性', '无', 1.0, N);

子图(1,3,3)；imagesc(真实(I3))； title('低分辨率过滤')
颜色图(灰色(256))；轴图像；轴关闭