Search code examples
matlabplotspectrumfrequency-analysis

Matlab power spectrum analysis


I would like to reproduce this image, but with my own EEG data. For what I understand, it is a power spectrum analysis done on filtered data.

single-sided amplitude spectrum after filtering

I recorded the EEG signal with a sampling rate of 1000Hz, with DC amplifiers (Low: DC; High:200). My Data is: 68 (electrodes) x 185080 (data points).

I have tried to use the following code from:http://uk.mathworks.com/help/signal/ug/psd-estimate-using-fft.html

 Fs = 1000;
 t = 0:1/Fs:1-1/Fs;
 x = Data;
 %x = cos(2*pi*100*t) + randn(size(t));

 N = length(x);
 xdft = fft(x);
 xdft = xdft(1:N/2+1);
 psdx = (1/(Fs*N)) * abs(xdft).^2;
 psdx(2:end-1) = 2*psdx(2:end-1);
 freq = 0:Fs/length(x):Fs/2;
 plot(freq,10*log10(psdx))
 grid on
 title('Periodogram Using FFT')
 xlabel('Frequency (Hz)')
 ylabel('Power/Frequency (dB/Hz)')

but this is what I obtain:

Periodogram using FFT

I am struggling understanding how to proceed in order to obtain an analysis of my EEG signal as in the first image. Any help would be very much appreciated.


Solution

  • Here's a simple example of how to do PSD using fft from scratch and without the DSP toolbox:

    %this does not include any filtering
    
    x = [0:0.01:pi];
    y = sin(100*x);
    nfft = 2^nextpow2(length(y));
    Fs = 100;
    
    psd1 = abs(fft(y,nfft)).^2/length(y)/Fs;%compute the PSD and normalize
    
    plot([0:50/(length(psd1)/2):50],psd1(1:length(psd1)/2+1))
    xlabel('Frequency (Hz)');
    ylabel('PSD');
    grid on
    title('PSD from FFT');
    

    The result:

    enter image description here

    If this method results in a similar one to what you posted then I think the other's comments about your data having some issues may be valid.