Possible fine-detail error in Matlab fft?

5 次查看（过去 30 天）
Tran Trung Luu
Tran Trung Luu 2021-3-15
评论： David Freiman 2022-2-19
Dear all,
I often work with FFT and I like the fft function from Matlab a lot:
https://www.mathworks.com/help/matlab/ref/fft.html
However, recently I encoutered a possible problem with the current implementation of Matlab:
In the standard calculation of the frequency axis that I extracted from the above link:
f = Fs*(0:(L/2))/L;
Fs is the sampling frequency and L is the length of the signal, it seems that the frequency step (the distance between two consecutive frequencies in the frequency axis) is Fs/L. This is usually not a problem if L is large enough, for example in this case. However, when I deal with small L, for example L = 5, then I realize that the frequency step is actually not the inverse of the total time window, which it should theoretically be. If I want to do this theoretically perfect, I had to use Fs/(L-1) but then this is not compatible with how the fft is done in Matlab, i.e. I got inconsistency of frequency peaks' position. It seems to me that the fft was not implemented with "perfect theoretical consideration".
Could you please advise/comment? Or you think I missed something?
Thanks and best regards,
Trung

请先登录，再回答此问题。

回答（1 个）

Prabhan Purwar
Prabhan Purwar 2021-3-21
Hi,
As the frequency resolution is defined as Fs/L, both the implementation and example seem rightly scaled and working as expected. If you can provide the calculations and steps you have tried, documents you have referred it will be helpful in resolving and understanding the issue more clearly.
Thanks
  2 个评论
Tran Trung Luu
Tran Trung Luu 2021-3-25
Hi Prabhan,
Thanks for your message.
I fully understand that "both the implementation and example seem rightly scaled and working as expected". I said this is a potential "fine-detail" problem. The problem is not big, but the implementation is probably not perfect, and that's why I want to discuss to find the perfect solution.
Here is a more detailed description:
Theoretically speaking, fft in Matlab is the fast DFT implementation of the Fourier Transform. From the intrinsic formalism for Fourier Transform, we know that the frequency step (df) is reciprocal of the time window (T) and vise versa, the frequency window (F) is the reciprocal of the time step (dt). If you consider now a discretized time axis or frequency axis, they are reprensented by N points, from t_1 to t_N, and the time step is dt = t_2 - t_1. Similarly we have f_1 to f_N and df = f_2 - f_1;
The standard description of frequency step/resolution in Matlab is Fs/L. This is not really wrong, but let's consider the following case: if we have only three points, N = 3. In this case, if you take L = N = 3, then the calculated fft result matches well the frequency grid calculated using Fs/L. However, it is physically wrong. The reason is that, in this case, the time window is (N-1)*dt which is only 2dt, not 3dt. So the actual frequency step, which it should theoretically be, should be 1/(2dt). But the real one calculated with the given formula is 1/(3dt).
In short, the discrepancy comes from the small number of the samples. If L is more than 1000, for example, then it does not matter much if we do Fs/1000 or Fs/999. For engineering purposes, it is ok in most cases. However, for scientists like me, it is important to make sure it is perfect.
Let me know how you think.
Thanks and best regards,
Trung
David Freiman
David Freiman 2022-2-19
Hi Trung,
I recently had a similar question and came across this post. I believe the resolution to this issue lies in the fact that for a discrete time signal as you describe above, with N=3, the time window is actually 3 * dt, because each sample point represents a time of dt.
In the past I have made the mistake of calculating the time window as t(N) - t(1). As you pointed out above, this would barely be noticeable for a large N of, say, 1000. However, when recently checking some code with a very low number of samples, N, I realized that the time window should be calculated as:
T = N * dt = N / Fs
So delta F would be 1 / T = Fs / N
Please let me know if you have any questions or see any errors in my logic.
Thank you,
Dave

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Transforms 的更多信息

标签

产品


版本

R2019b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by