clear all;
dir='/home/azim/mlk/data/meg/exp1a/';

cd(dir)

subj = {
'fiona/LH_IFG_VE_fi'
'hannah/LH_IFG_VE_ha'
'ian/LH_IFG_VE_ih'
'ianf/LH_IFG_VE_if'
'james/LH_IFG_VE_ja'
'jenni/LH_IFG_VE_je'
'michael/LH_IFG_VE_mi'
'ned/LH_IFG_VE_ne'
'pete/LH_IFG_VE_pe'
};

N=max(size(subj));
freq=cell(1,N);

% INDUCED PLUS EVOKED %
% PRE-ROCESSING PER SUBJECT %

for i=1:N

   sprintf('SUBJECT %d\n',i)
   dataset = [dir subj{i} '.ds'];

   cfg1 = [];
   cfg1.dataset = dataset;
   cfg1.channel = {'MEG' 'STIM' 'V0'};
   cfg1.trialdef.excludeConditions = {'BAD'};
   cfg1.trialdef.eventtype = 'tr15';
   cfg1.trialdef.prestim   = 0.25;
   cfg1.trialdef.poststim  = 0.7;
   cfg1.detrend            = 'yes';
   cfg1.blc                = 'yes';
   %cfg1.keeptrials        = 'yes';

   cfg1 = definetrial(cfg1);
   raw{i} = preprocessing(cfg1);

end

% INDUCED PLUS EVOKED %
% FREQUANALYSIS PER SUBJECT - KEEPING TRIALS %

for i=1:N

   cfg2 = [];
   cfg2.output       = 'pow';
   cfg2.sgn          = 'MEG';	
   cfg2.method       = 'wltconvol';
   cfg2.channel      = 'V0';
   cfg2.keeptrials   = 'yes';
   cfg2.width        = 3;                     % Width of wavelet
   cfg2.gwidth       = 2.5;
   cfg2.foi          = 1:0.2:40;	
   cfg2.toi          = -0.25:0.005:0.7;
   iefreq{i} = freqanalysis(cfg2, raw{i});

end

% INDUCED PLUS EVOKED %
% COMPUTE AVERAGE POWER FOR EACH SUBJECT %

for i=1:N

  cfg3 = [];
  cfg3.foilim        ='all';
  cfg3.toilim        ='all';
  cfg3.channel       ='V0';
  miefreq{i} = freqdescriptives(cfg3, iefreq{i});

end

% INDUCED PLUS EVOKED %
% BASELINE NORMALISATION FOR T_F DATA PER SUBJECT %

for i=1:N

   cfg4 = [];
   cfg4.baseline      =[-0.25 0.0];
   cfg4.baselinetype  ='absolute';
   nmiefreq{i} = freqbaseline(cfg4,miefreq{i});

end

% INDUCED PLUS EVOKED %
% PLOTTING INDIVIDUAL SPECTROGRAMS %

cfg5 = [];
cfg5.baseline     = 'no';
%cfg5.baseline     = [-0.25 0.0];  	
%cfg5.baselinetype = 'relative';
%cfg5.baselinetype = 'absolute';
cfg5.zlim         = 'maxmin';
cfg5.xlim         = [0.0 0.7];
cfg5.ylim         = [1 45]; 	
cfg5.channel      = 'V0';

for i=1:N
  figure(22);singleplotTFR(cfg5, nmiefreq{i});
  pause
end

% INDUCED PLUS EVOKED %
% COMPUTING FREQ GRANDAVERAGE %

cfg6=[];
cfg6.keepindividual = 'yes';
cfg6.foilim         = 'all';
cfg6.toilim         = 'all';
cfg6.channel        = 'V0';

ind_iefgrandavg=freqgrandaverage(cfg6, nmiefreq{1}, nmiefreq{2}, nmiefreq{3}, nmiefreq{4}, nmiefreq{5}, nmiefreq{6}, nmiefreq{7}, nmiefreq{8}, nmiefreq{9});



% EVOKED %
% TIMELOCKED AVERGAING PER SUBJECT %

for i = 1:N;

   cfg8 = [];
   cfg8.channel       = 'V0';
   cfg8.latency       = [-0.25 0.7];
   cfg8.keeptrials    = 'no';
   cfg8.normalizevar  = 'N-1';
   tav{i} = timelockanalysis(cfg8,raw{i});

end

% EVOKED %
% FREQUENCY ANALYSIS  PER SUBJECT %

for i=1:N

   cfg9 = [];
   cfg9.output       = 'pow';
   cfg9.sgn          = 'MEG';	
   cfg9.method       = 'wltconvol';
   cfg9.channel      = 'V0';
   cfg9.keeptrials   = 'no';
   cfg9.width        = 3;                     % Width of wavelet
   cfg9.gwidth       = 2.5;
   cfg9.foi          = 1:0.2:40;	
   cfg9.toi          = -0.25:0.005:0.7;
   efreq{i} = freqanalysis(cfg9, tav{i});

end

% EVOKED %
% BASELINE NORMALISATION FOR T_F DATA PER SUBJECT %

for i=1:N

   cfg10 = [];
   cfg10.baseline      =[-0.25 0.0];
   cfg10.baselinetype  ='absolute';
   nefreq{i} = freqbaseline(cfg10,efreq{i});

end


% EVOKED %
% PLOTTING INDIVIDUAL SPECTROGRAMS %

cfg11 = [];
cfg11.baseline     = 'no';
%cfg11.baseline     = [-0.25 0.0];  	
%cfg11.baselinetype = 'relative';
%cfg11.baselinetype = 'absolute';
cfg11.zlim         = 'maxmin';
cfg11.xlim         = [0.0 0.7];
cfg11.ylim         = [1 45]; 	
cfg11.channel      = 'V0';

for i=1:N
  figure(23);singleplotTFR(cfg11, nefreq{i});
  pause
end


% EVOKED %
% COMPUTING FREQ GRANDAVERAGE %

cfg12=[];
cfg12.keepindividual = 'yes';
cfg12.foilim         = 'all';
cfg12.toilim         = 'all';
cfg12.channel        = 'V0';

ind_efgrandavg=freqgrandaverage(cfg12, nefreq{1}, nefreq{2}, nefreq{3}, nefreq{4}, nefreq{5}, nefreq{6}, nefreq{7}, nefreq{8}, nefreq{9});



% EVOKED VERSUS EVOKED+INDUCED%
% STATISTICAL COMPARISON %


cfg15 = [];
cfg15.statistic = 'depsamplesT';
cfg15.alphathresh = 0.05;
cfg15.makeclusters = 'yes';
cfg15.minnbchan = 2;
cfg15.clusterteststat = 'maxsum';
cfg15.onetwo = 'twosided';
cfg15.alpha = 0.05;
cfg15.nranddraws = 100;
cfg15.channel = {'V0'}
%cfg15.latency = [0 0.5]  % increases sensitivity if prior temporal information is knowm
cfg15.frequency = [5 40];

%frange = 'all_freqs';  % file addition for saving

[clusrand] = clusterrandanalysis (cfg15, ind_efgrandavg, ind_iefgrandavg);