Using fft() to calculate PSD.

I used the code below perform spectral analysis of signal with known period of 4 s:

load lujan.dat
data_out.time = lujan(:,1);
data_out.disp = lujan(:,2);

Y = fft(data_out.disp,1000);
Y(1)=[];

NY = length(Y);
Pyy = abs(Y(1:floor(NY/2))).^2;
nyquist = 1/2;
f = (1:NY/2)/(NY/2)*nyquist;

plot(f,Pyy,'-b','LineWidth',3);
title('Power spectral density');
xlabel('Frequency (Hz)');

However the output graph shows a peak at 0.05 Hz which is not correct! IS there something obvious i'm doing wrong?

thanks

Answer 1

[Wayne King]

Are you sure the sampling frequency is 1 ? If so, how far off is the peak from where you expect it to be? Is it almost correct or way off? You also don't tell us how long the time series is, so we don't know why you are adding the length input argument to fft() . Your frequency increment of 1/NY (0.0025) is correct.

[2one]

the signal is from 0 to 12 sec, sampled every 0.2 sec (61 points total).

the signal has period of 4 s so frequency peak should be at 0.25 hZ but graph is showing peak at 0.05 Hz (which would mean 20 s period).

Answer 2

[Wayne King]

Then your sampling frequency is specified incorrectly.

 Fs = 1/0.2;
 % here I'll create a signal just to show you
 t = 0:1/Fs:12-1/Fs;
 x = cos(2*pi*0.25*t)+randn(size(t));
 Y = fft(x,256);
 Pyy  = abs(Y(1:floor(256/2)+1)).^2;
 freq = 0:Fs/256:5/2;
 plot(freq,10*log10(Pyy))
 grid on
 set(gca,'xtick',[0 0.25 0.5 1 1.5 2 2.5])