博客

用MATLAB进行图像处理

图像处理的概念、算法和MATLAB

FFT光谱泄漏

发布的史蒂夫·埃丁，

333 views (last 30 days 喜欢| 3评论
一个MATLAB用户最近联系了MathWorks技术支持，问为什么输出金宝app FFT. 没有达到他们的期望，技术支持要求Matlab数学团队获得帮助。金宝app乔治亚州科技毕业生克里斯转衣写了一个我喜欢阅读的详细回应。我认为值得适应博客帖子的情况。虽然这种情况是关于1-D FFT，但潜在的问题也出现在图像处理中，我已经在过去写了一下。（见我的 傅里叶变换的类别 ．）
让我们从一个50hz的正弦波开始，采样频率为1khz，采样频率为4080个样本。
FS = 1000;
N = 4080;
t = 0: (n - 1) / Fs;
y =罪(2 *π* 50 * t);
情节（t，y）
轴（[0 0.1 -1.2 1.2]）
Xlabel（“t (sec)”）
ylabel（'y（t）'）
网格在
标题（'y（t）'的第一个0.1秒）
接下来，让我们计算和显示FFT，缩放频率轴，使其处于Hz，并通过FFT长度的平方根缩放幅度。
Y = fft (Y);
f =（fs / n）*（0：n-1）;
绘图（F，ABS（Y）/ SQRT（n））
Xlabel（的频率(赫兹)）
网格在
让我们在50赫兹围绕峰值放大。
XLIM（[45 55]）
到目前为止，这似乎是合理的。现在，让我们只使用前4096个样本尝试该过程 $ y（t）$ ．（大约4.1秒。）
N2 = 4096;
t2 = (0: (N2-1)) / Fs;
y2 = sin（2 * pi * 50 * t2）;
Y2 = fft (Y2);
f2 = (Fs/N2) * (0:N2-1);
情节(f2、abs (Y2) /√(N2))
Xlabel（的频率(赫兹)）
XLIM（[45 55]）
网格在
看起来很多。让我们把它们放在一起。
持有在
绘图（F，ABS（Y）/ SQRT（n））
持有从
传奇([“| y_2 |”，“Y | |”]）
嗯。带有4080个点的信号的FFT在50hz处有一个清晰的、尖锐的峰值，而带有4092个点的信号的FFT在那个峰值附近分布得更广。
这是技术支持问题：差异的解释是什么？金宝app有什么问题 FFT. 功能 $ n = 4092 $ ？
克里斯立即给出了一个很好的，概念性的解释，我马上就会讲到。不过，首先，我想回顾一下我在遇到类似问题时有时会用到的思考过程。是否有一种独立的方法来验证问题中的函数是否返回了正确的答案?在这种情况下，我知道 FFT. 功能计算一个名为the的东西 离散傅里叶变换 或DFT。对于离散时间信号 $ y [n] $ ，定义为 $ 0 \leq n < n $ ， DFT为:
$ Y [k] = \ sum_ {n = 0} ^ {n} Y [n] e ^ {- j 2 \πk n / n}, k = 0, 1 \ ldots, n - 1美元
的 FFT. 函数不使用此公式计算DFT，因为它很慢，但我们可以使用公式来仔细检查其输出。
k = (0: N2-1);
n = k';
y2p = sum（y2'。* exp（-1i * 2 * pi * k。* n / n2））;
Y2. 和 Y2p 与相对精度的约12个小数位内相同：
规范(Y2p - Y2) /规范(Y2)
ans = 7.9563e-13
这在普通浮点圆差异中，我期望看到的差异很大。但是让他们绘制它们，只是为了确定。
情节(f2、abs (Y2) /√(N2), f2, abs (Y2p) /√(N2))
网格在
XLIM（[45 55]）
这两个结果在视觉上难以区分。
好的，此时，那么，我会相信这一点 FFT. 返回正确的答案。我们可以转向为什么输出似乎如此不同的问题，因为我们稍微改变了样本总数。让我把它交给克里斯，谁提供了以下解释：
这只是 光谱泄漏 ;这是众所周知的 很好的描述 傅里叶分析中的问题。由长度的DFT采样的离散光谱频率 N 是 (0: (N - 1) / N ．在客户的情况下，他的纯音的离散频率不是整数倍数 $ 1 / n $ ．因此，频率内容“泄漏”到其他容器中。
看到这一点的另一种方式是回顾说明DFT / FFT固有地假设输入是周期性信号。但是，既不上面的信号， y 和 Y2. ，包含整数个数的句点。你可以看到如果你画的是信号的结束而不是开始。
情节(t ((end-50):结束),y ((end-50):结束)
ylim（[ - 1.2 1.2]）
标题（“N = 4080”）
网格在
情节(t2 ((end-50):结束),y2 ((end-50):结束)
ylim（[ - 1.2 1.2]）
标题（“N = 4096”）
网格在
与 $ n = 4080 $ ，信号在结束时几乎完成了它的周期，所以最后一个样本几乎等于第一个样本(它是0)。因此，假设的周期性产生的伪信号并不强。虽然有4096个样本，但信号结束时仍有四分之一的周期。从上一个样本值(-1)到第一个样本值(0)有一个大的跳跃，这个大的跳跃就像信号中的一个不连续，使频谱峰扩散，导致更多的“泄漏”。
如果一切顺利，那么生活变得美丽;与任何整数 $ k \ geq 2 $ ， 使用 $ n = 50k $ 会产生我们都想看到的理想峰值。然而，在其他方面 $ 50/1000美元 不能写成 $ q / n $ 为一个整数 问 ，就会出现频谱泄漏现象。振幅会下降，峰值会扩散。结果是正确的;这是数学上的预期结果，即使它不符合我们最初的直觉。
谢谢你让我用你的解释，克里斯!
|
类别：
傅里叶变换
  • 打印
  • 发送电子邮件

另请参阅

注释

要发表评论，请点击这里登录您的MathWorks帐户或创建新的。