MATLAB: Discrepancies between Power Spectra and Data

psdpwelch

Hi MATLAB users,

I have time-series of repetitive burst, like the following link (lower trace): http://mitpdev.mit.edu/library/erefs/arbib/images/figures/A002_fig001.gif

The bursting occurs at a frequency of 0.5Hz, and within the burst there are frequencies of between 50-70Hz. Sampling at 10millisecs, my power spectra;

Fs=100;
pwelch(event_train,128,120,[],Fs,'onesided')

doesn't give me what I would expect (peaks at 0.5Hz and 50-70Hz). What am I doing wrong?

Thanks,

Linford

Best Answer

Hi Linford, A sampling frequency of 100 Hz (dt=0.01) is going to result in aliasing for the content from 50-70 Hz. As far as your other question, is there anyway you can post your time series anywhere? How long is your data vector? Your frequency resolution with a segment length of 128 points and your sampling frequency is almost 0.4 Hz so one thing to do would be to increase your segment length.

Without seeing the data or your PSD estimate, it's hard to say anything more. Another thing is what is the mean of the data. Is the mean nonzero?

If so, the power at DC (zero frequency) can obscure the spectral details you are interested in, so you may want to detrend() first.

Related Solutions

MATLAB: Autocorrelations and Correlogram

Hi Linford, If you use xcorr() with the 'coeff' option you get the normalized autocorrelation sequence which does not exceed 1 in absolute value.

The triangle shape that you speak about is expected. The sample autocorrelation sequence by necessity uses fewer and fewer samples at larger lags (the summation cannot exceed the length of the input vector), so you do expect the sample ACF to decay.

Here is an example with a sinusoid. I'll use a 100 Hz sine wave sample at 1 kHz.

t = 0:0.001:1-0.001;
x = cos(2*pi*100*t);
[xc,lags] = xcorr(x,'coeff');
stem(lags,xc)

Now let's zoom in a bit.

stem(lags(length(x):length(x)+50),xc(length(x):length(x)+50));
xlabel('Lags'); ylabel('ACF');

Now you see that the ACF is oscillating with a period of 10 lags, which is exactly what you expect since the sampling interval is 1 msec and the period is 10 msec.

MATLAB: Looking for low frequncy peaks from the FFT of a time series signal

Your frequency peaks are much smaller than the d-c offset in ‘dsignal’, so the d-c offset will dominate the Fourier spectrum. The easiest way to eliminate it is to subtract the d-c offset (mean value of your signal) from the signal before taking the Fourier transform of it.

Try this:

L = numel(dsignal);
Fs = 25;                                                % Sampling Frequency (Hz)
Fn = Fs/2;                                              % Nyquist Frequency (Hz)
t = linspace(0, 1, L)/Fs;                               % Time Vector
y = fft(dsignal - mean(dsignal))/L;                     % Fourier Transform, Normalise For Length
Fv = linspace(0, 1, fix(L/2)+1)*Fn;                     % Frequency Vector (One-Sided Fourier Transform)
Iv = 1:numel(Fv);                                       % In \dex Vector
figure
plot(Fv, abs(y(Iv))*2)
grid
xlim([0 0.5])

That will clearly show the 3 low-frequency peaks, and some of the low-frequency noise in the rest of the signal.

Best Answer

Related Solutions

MATLAB: Autocorrelations and Correlogram

MATLAB: Looking for low frequncy peaks from the FFT of a time series signal

Related Question