-
Notifications
You must be signed in to change notification settings - Fork 7
/
doFlynn.m
49 lines (36 loc) · 2.1 KB
/
doFlynn.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
function [fftResults,frequencies] = doFFT(trimmedEEG)
% doFFT run a fft on EEG data
% based on code by Michael X. Cohen
% modified by C. Hassall and O. Krigolson
% May 23, 2018 last update
% The input data is channels x time x epochs (3 D matrix)
% The output data is in power (amplitude squared)
% create a Hanning Window to reduce edge artifacts - divide by two to get
% half a cosine cycle - the discontinuity from edge artifacts is
% because the start and end voltages of the segment are different
% Mike Cohen says to turn off the Hanning Window
% hanning_window = [];
% hanning_window = 0.5 - cos(2*pi*linspace(0,1,length_of_data))/2;
% convolute the data with the hanning window
% tapered_data = trimmedEEG.data .* hanning_window;
% run the fft
fftOutput = fft(trimmedEEG.data,[],2);
% determine frequency resolution
frequencyResolution = trimmedEEG.srate / trimmedEEG.pnts;
%We want to scale the output of the fft function to be divided by the length of the data submitted because The FFT computation involves a lot of summing over time points, and so dividing by N puts the data back to the scale of the original input.
fftOutput = fftOutput/trimmedEEG.pnts;
% we want to remove the negative half of the fft output because it is redundant
fftOutput = fftOutput(:,1:trimmedEEG.pnts/2,:);
% you take abs of the output of the fft because you want the resultant of the real and imaginary output
fftOutput = abs(fftOutput);
% square the amplitude to get power
fftOutput = fftOutput.^2;
% multiply by 2 to ensure that the positive aspect of the fft output has been corrected for the loss of the negative aspect
fftOutput = 2 * fftOutput;
% finally, remove the very first value as it is the power at 0 Hz (the DC power)
fftOutput(:,1,:) = [];
% Compute the average FFT output for this condition
fftResults = squeeze(mean(fftOutput,3));
% create a frequency bin
frequencies = frequencyResolution:frequencyResolution:length(fftOutput)*frequencyResolution;
end