3 views (last 30 days)

Show older comments

Hello.

In the Matlab description of the FFT implementation (https://www.mathworks.com/help/matlab/ref/fft.html). There is this part where the vector Y is divided by the length L in order to compute the spectrum. I cannot understand why it is required to divide by L, if someone can maybe enlighten this point.

Compute the two-sided spectrum P2. Then compute the single-sided spectrum P1 based on P2 and the even-valued signal length L.

P2 = abs(Y/L);

P1 = P2(1:L/2+1);

P1(2:end-1) = 2*P1(2:end-1);

Thank you.

dpb
on 8 Jul 2021

In short, it's simply to normalize the PSD peak to the amplitude of the input time series -- if don't then the magnitude is proportional to the length of the input.

Just had long discussion (amongst many over the years) a short time ago. The final comment I added at https://www.mathworks.com/matlabcentral/answers/847665-why-do-i-receive-results-mismatch-in-the-fft-signal-while-using-matlab-built-in-fft-function-how-ca?s_tid=srchtitle#comment_1593465 illustrates what happens if one doesn't normalize;

David Goodmanson
on 9 Jul 2021

Hi Anas,

suppose you have N = 1000 points over an interval of one second. The time array is t = (0:999)/1000. [see note below]. Let f0 = 1 Hz and consider the complex waves

g_n(t) = exp(2*pi*i*n*f0*t) for n = 1,2,3, ....

Each g_n(t) oscillates n times in the time window.

The fft takes the signal, and for each m = 1,2,3 ... [1] multiplies the signal by g_m(t)* = exp(-2*pi*i*m*f0*t), [2] does a sum over all the array points and [3] reports out the answer.

Suppose your signal is a single oscillatory wave of amplitude 1 for a particular n0. After step [1] the fft will do a sum over all the points of the expression

B = exp(2*pi*i*(n0-m)*f0*t).

The result is 0 except when m = n0. In that case B = contant = 1 and the sum over the N points gives the result N. So if you want to recover the original amplitude 1 you have to divide the fft result by N.

[note] The array has 1000 points and 1000 intervals, including the interval from .999 sec to 1 sec, but does not include the point at 1 sec.

David Goodmanson
on 10 Jul 2021

Every signal can be represented as the sum of complex waves with amplitude coefficients cn:

sig = Sum{n} cn*exp(2*pi*i*n*f0*t)

The fft determines the cn coefficients, which is the frequency spectrum. To do that multiply the signal by g_m(t)* = exp(-2*pi*i*m*f0*t), take the sum over points and divide by N, which gives, on the right hand side,

Sum{n} cn*exp(2*pi*i*(n-m)f0*t) /N

AsI mentioned before, the sum is N when n = m and zero otherwise. You get N*cm/N = cm. Do that for each m. That's the principle. As for the underlying software, that's proprietory with Matlab, but if you look up the Cooley-Tukey algorithm you will get the basic idea.

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

Start Hunting!