Building an FFT on an octave

Question

I know that FFT changes the function in the time domain to one shown in the frequency domain.

However, when I try to plot in the frequency domain, I can get it to work normally, using time as the X axis, when it was supposed to be not the same, but the frequency.

In addition, I can only get the amplitudes so that they match the values ​​in the original signal by dividing the y axis by some integer. Why is this?

Here is my code

t=0:0.001:2 x=2*sin(20*pi*t) + sin(100*pi*t) subplot(2,1,1) plot(1000*t,x) grid xlabel("Time in milliseconds") ylabel("Signal amplitude") subplot(2,1,2) y=fft(x) plot(1000*t,abs(y)) xlabel("Frequency") ylabel("Signal amplitude") 

and graphics.

Please help = (

Answer 1

Frequency ratio (x-axis scaling)

The frequency of each value obtained by the FFT is linearly related to the index of the output value through:

 f(i) = (i-1)*sampling_frequency/N 

Where N is the number of FFT points (i.e., N=length(y) ). In your case, N=2001 .

You can subtract the sampling rate from your definition of t as 1 / T, where T is the sampling time interval (T = 0.001 in your case). Thus, the sampling rate is 1000 Hz.

Note that since the value of t(i) also linearly related to the index i , through

 t(i) = (i-1)*0.001 

it is possible (although it is not necessary to advise, as this will simply obscure your code) to determine f = 1000*t*sampling_frequency/N Please note that you are missing the sampling_frequency/N member, which accordingly led to the display of tones with the wrong frequency (from the definition x should be peaks with a frequency of 10 Hz and 50 Hz and corresponding aliases at frequencies of 990 Hz and 950 Hz).

Amplitude Ratio (Y-axis Scaling)

Please note that the observed ratio is only approximate, therefore the following is not a mathematical proof, but simply an intuitive way to visualize the relationship between tone amplitudes in the time domain and peak values ​​of the frequency domain.

Simplification of a problem to one tone:

 x = A*sin(2*pi*f*t) 

The approximate amplitude of the corresponding peak can be obtained using the Parseval theorem :

In the time domain (left side of the equation), the expression is approximately 0.5*N*(A^2) .

In the frequency domain (right side of the equation), making the following assumptions:

• spectral leakage effects are negligible.
• the spectral content of the tone is contained in only 2 pins (at the frequency f and the corresponding smoothing frequency sampling_frequency-f ) summation is taken into account (all other bins are equal to ~ 0). Note that this is usually only done if the tone frequency is an exact (or nearly accurate) multiple of sampling_frequency/N

the expression on the right-hand side is approximately 2*(1/N)*abs(X(k))^2 for some value of k corresponding to the peak at frequency f .

Putting the two together, we get abs(X(k)) ~ 0.5*A*N In other words, the output amplitude shows a scaling factor of 0.5*N (or approximately 1000 in your case) relative to the amplitude in the time domain, as you noticed.

The idea is still applied with more than one tone (although the negligible assumption of spectral leakage ultimately breaks down).

Answer 2

Other answers suggested that in this example there are frequency responses at 950 Hz and 990 Hz. This is a misunderstanding about how FFT code uses indexes. These "high frequency" peaks are actually -50Hz and -10Hz.

The frequency domain extends from -N / 2 * sampling_frequency / N to + N / 2 * sampling_frequency / N. But for historical reasons, it is accepted that the first N / 2 pieces of information are positive frequencies, the middle point is zero frequency, and the last N / 2 pieces of information are negative frequencies in reverse order. For the power spectrum, there is no need to show more than the first 1 + N / 2 pieces of information.

This agreement is extremely confusing since I had to puzzle it from Press et al. Hartley Fast Transform numerical recipes and manual coding, many years ago when I first used FFT, even before the Matlab 1.0 beta test, which Cleve Mohler handed out to some lucky ones :-)

Answer 3

You can use fftshift to center at zero frequency, but this is not so important.

According to the Nyquist / Kotelnikov theorem, the highest frequency you get is half the sampling frequency: fmax = fsample / 2. Since you have a sampling rate of 1 kHz, you get a maximum of 500 Hz. The frequency step is the inverse of the reception time: df = 1/2 = 0.5 Hz in your case.

So, the formula indicated in the other answer is incorrect, you should double it.

You now have 1000 in the middle of the graph. You will use only the first half. For each point, you get half the number on the x axis to get the frequency.

PS.I could not comment, because I do not have enough reputation.