Power spectra

A couple of places calculate a power spectrum. The process is as follows.

For the input timeseries y:

Demean

 ynorm = y-mean(y);

Calculate the power spectral density with “periodogram”

fs = 1/TR;
[pxx,f] = periodogram(ynorm,[],[],fs);

Remove meaningless low frequencies

% Remove meaningless low frequencies.
minf = 1/(t(end)/2);
indCutoff = find( f < minf , 1, 'last' );

pxx(1:indCutoff) = [];
f(1:indCutoff) = [];

Calculate total power

    powerTotal = bandpower( pxx, f, 'psd' );