MATLAB: DFT using matlab function

digital image processingfftfrequencysignal processing

In case a off-nominal signal is processed using nominal frequency based fft algorithm, the peak value of the magnitude of fft is still at the frequency index or bin pertaining to the nominal value (though the amplitude is different and the side lobes also have an amplitude which is not seen when the signal has a nominal value of frequency). Why does the frequency index not change even though the frequency is no longer the nominal value. How can i visualize this in matlab using the fft function.

Best Answer

Hi mk,

"Why does the frequency index not change even though the frequency is no longer the nominal value "

When you do an fft, the frequency array is predetermined by the properties of the time array. For an N-point fft and array spacings of delt and delf, the golden rule for ffts is

delt*delf = 1/N.

regardless of the properties of the signal. In the following code, delf = 1 so the frequency array represents exact integral frequencies. The example uses a cosine wave with f0 = 10. The total time record is 1 sec and you can see exactly 10 cycles in the time plot. In the frequency plot, fftshift is used to put f = 0 at the center of the plot. There are two peaks exactly at +-10Hz with amplitude 1/2, illustrating that

cos(2*pi*f0*t) = ( exp(i*2*pi*f0*t) + exp(-i*2*pi*f0*t) )/2.

and no frequency content anywhere else.

Now change f0 to 10.2 Hz and run it again. The frequency grid is incapable of exactly representing 10.2 Hz so it gives a reduced peak at 10 Hz, with frequency content spilling over into other frequencies.

N = 1e4;
delt = 1e-4         % total time record = 1 sec
delf = 1/(N*delt);  % delf = 1;
          % or you could do what is commonly done
          % Fs = 1e4; sampling rate
          % delt = 1/Fs;
          % delf = Fs/N;
t = (0:N-1)*delt;
f = (-N/2:N/2-1)*delf;    % freq grid for fftshift
f0 = 10;
y = cos(2*pi*t*f0);
figure(1)
plot(t,y); grid on
z = fftshift(fft(y)/N);   % ordinarily z is complex, real in this case
figure(2)
plot(f,z,'o-'); grid on
xlim([-40,40])

Related Solutions

MATLAB: I wish to get the spectrum (Magnitude vs Frequency(Hz)) of a pure tone(1 kHz). However, the piece of code I’m employing to do so isn’t working. What am I doing wrong

Hello Komal,

The whole point of using fftshift is to put zero frequency at the center of the resulting frequency array. If you zero-pad fft_puretone on one side only, it messes up that symmetry. So fft_puretone needs to be zero-padded equally on each side.

An an N-point fft with time array spacing of delt and frequency array spacing of delf always obeys delt*delf = 1/N. That remains true after subsequent zero-padding.

For the correct frequency domain amplitude, you need to multiply the fft by a factor of 1/N.

You did not specify what puretone looks like in the time domain, but I am speculating that the time record was 5 seconds, and the amplitude of the wave = 1/sqrt(2) (peak), = 1/2 (rms).

Here is an example of how this would work. The resulting arrays have 1105920 elements and show a peak at 1kHz and another at -1kHz. After ffting a real function, people often double the positive frequency peak and ignore the negative frequency peak, which is all right as long as you keep in mind what the fft is actually doing.

hmax = 1105920
fs = 44100
% 1k puretone in time domain
  delt = 1/fs;
  f0 = 1000;
  T = 5;             % time record (sec)               
  N = T*44100;
  t = (0:N-1)*delt;
  puretone = (1/sqrt(2))*sin(2*pi*f0*t);  
% back to fft code
  fft_puretone = abs(fftshift(fft(puretone)/N));
  Npad = hmax-N;
  fft_puretone = [zeros(1,Npad/2) fft_puretone zeros(1,Npad/2)];
  delf = 1/(N*delt);   % fundamental rule for fft
  freq = (-hmax/2:hmax/2-1)*delf;
  figure(1)
  plot(freq,fft_puretone)

MATLAB: FFT and PSD – normalize values

Hi John,

This concerns fft only; psd is a different animal.

For a continuous oscillation, most people do not seem to be concerned about fft output that spills into adjacent frequencies, but it's worth looking at as you are doing. With the fft you can get sharp one-point frequency peaks with a varying number of points and a varying amplitude, but not for time windows of arbitrary width. With fft It's required that there be an exact integral number of oscillations in the time window. And, using a single oscillation of a sine wave as an example, the last point in the array can't repeat the first point. The last point has to be "one point short" of a full oscillation

y = sin(2*pi*(0:9)/10)    % correct   but   y = sin(2*pi*(0:10)/10)      % incorrect

The code below gives two peaks, one at +160 kHz and one at -160 kHz for any number of fft points (above a minimum number), any number of complete oscillations in the time window and any amplitude. The abs(amplitude) of each peak is A/2, consistent with cos(w*t) = (1/2)(exp(i*w*t)+exp(-i*w*t)). If you want, as is common for plotting purposes, you could ignore the negative freqencies and double the height of the positive frequency peak.

f0 = 160e3;
ncyc = 23;     % number of cycles in time window
A = 6;         % amplitude
psi = pi/6;    % arbitrary phase angle
N = 1299;      % need N >2*ncyc to avoid aliasing
t0 = 1/f0;
T = ncyc*t0;   % time window
delt = T/N;
t = (0:N-1)*delt;
y = A*cos(2*pi*f0*t+psi);
figure(1)
plot(t,y);grid on
delf = 1/(delt*N);    % golden rule for fft
if rem(N,2) ==0
  f = (-N/2:N/2-1)*delf;
else
  f = (-(N-1)/2:(N-1)/2)*delf;   
end
z = fftshift(fft(y))/N;     % fftshift to put f=0 in the middle
figure(2)
stem(f,abs(z),'o-');grid on

Best Answer

Related Solutions

MATLAB: I wish to get the spectrum (Magnitude vs Frequency(Hz)) of a pure tone(1 kHz). However, the piece of code I’m employing to do so isn’t working. What am I doing wrong

MATLAB: FFT and PSD – normalize values

Related Question