%% Fieldtrip TFA
clear all
close all

addpath('/Users/patrickwiegel/Documents/MATLAB/Scripts/2020/TF_Data/fieldtrip-20200911')
addpath('/Users/patrickwiegel/Documents/MATLAB/Scripts/2020/TF_Data/fieldtrip-20200911/external/eeglab')

cd '/Users/patrickwiegel/Documents/MATLAB/Scripts/2020/TF_Data/TF_Data'


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Prepare data for TF analyses


i=1;
load(['Feedback_raw_' num2str(i) '.mat'])
load(['EEG_' num2str(i) '.mat'])

load('Index_correct.mat')
load('Index_incorrect.mat')


%% import EEG data
EEG_Data = eeglab2fieldtrip( EEG, 'raw', 'none' );

%% define trials for frequency analyses for contrast 
cfg = [];
cfg.trials = Indexes_correct{1, 1};
cfg.toilim = [0.25 0.55];
data_fed_correct = ft_redefinetrial(cfg, EEG_Data);


cfg = [];
cfg.trials = Indexes_incorrect{1, 1};
cfg.toilim = [0.25 0.55];
data_fed_incorrect = ft_redefinetrial(cfg, EEG_Data);

%% run frequenc analyses for two contrasting conditions 
cfg = [];
cfg.method    = 'mtmfft';
cfg.output    = 'powandcsd';
cfg.tapsmofrq = 4;
cfg.foilim    = [25 25];
frequ_correct = ft_freqanalysis(cfg, data_fed_correct);

cfg = [];
cfg.method    = 'mtmfft';
cfg.output    = 'powandcsd';
cfg.tapsmofrq = 4;
cfg.foilim    = [25 25];
frequ_incorrect = ft_freqanalysis(cfg, data_fed_incorrect);

cfg = [];
cfg.method    = 'mtmfft';
cfg.output    = 'powandcsd';
cfg.tapsmofrq = 4;
cfg.foilim    = [25 25];
frequ_all = ft_freqanalysis(cfg,EEG_Data);





%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Source analyses starts 

%% load MRI and BEM templates
load('/Users/patrickwiegel/Documents/MATLAB/Scripts/2020/TF_Data/fieldtrip-20200911/template/headmodel/standard_bem.mat')
load('/Users/patrickwiegel/Documents/MATLAB/Scripts/2020/TF_Data/fieldtrip-20200911/template/headmodel/standard_mri.mat')

%% 
cfg = [];
mri = ft_volumereslice(cfg, mri);
display(mri)

%% convert coordinate system

cfg = [];
cfg.method = 'interactive';
cfg.coordsys = 'ctf';
mri_realigned = ft_volumerealign(cfg,mri);
display(mri_realigned)

%%
cfg = [];
cfg.method = 'flip';
mri_realigned = ft_volumereslice(cfg,mri_realigned);

%% segment mri 

cfg = [];
cfg.write      = 'no';
[segmentedmri] = ft_volumesegment(cfg, mri_realigned);

disp(segmentedmri)

%% align electrodes
load('elec_aligned')
elec_aligned = ft_convert_units(elec_aligned, 'cm');
vol = ft_convert_units(vol, 'cm'); % Convert the vol to cm, because the CTF convenction is to express everything in cm.

%%
[headmodel, elec_aligned]=ft_prepare_vol_sens (vol,elec_aligned);

cfg           = [];
cfg.method    = 'interactive';
cfg.elec      = elec_aligned;
cfg.headshape = headmodel.bnd(1);
elec_aligned  = ft_electroderealign(cfg);

%%
% construct the dipole grid in the template brain coordinates
% the negative inwardshift means an outward shift of the brain surface for inside/outside detection
cfg              = [];
cfg.resolution   = 1;
cfg.tight        = 'yes';
cfg.inwardshift  = -1.5;
cfg.headmodel    = headmodel;
sourcemodel    = ft_prepare_sourcemodel(cfg);

%%
cfg              = [];
cfg.method       = 'dics';
cfg.frequency    = 25;
cfg.sourcemodel  = sourcemodel;
cfg.headmodel    = headmodel;
cfg.elec         = elec_aligned;
cfg.dics.projectnoise = 'yes';
cfg.dics.lambda       = '5%';
cfg.dics.keepfilter   = 'yes';
cfg.dics.realfilter   = 'yes';
sourceAll = ft_sourceanalysis(cfg, frequ_all);





%%
%previous try

%% 
% cfg                  = [];
% cfg.headmodel        = headmodel;
% cfg.elec             = elec_aligned;
% cfg.sourcemodel  = sourcemodel;
% cfg.reducerank       = 2;
% cfg.resolution       = 1;   
% cfg.sourcemodel.unit = 'cm';
% sourcemodel = ft_prepare_leadfield(cfg,frequ_all);