From e948ab44798cdc9ef63ddeaf2e71a283c7a0b3de Mon Sep 17 00:00:00 2001
From: Dung Truong <dt.young112@gmail.com>
Date: Tue, 3 Oct 2023 16:31:15 -0700
Subject: [PATCH] remove avg reref from preprocessing but keep in plotting

---
 eeg_nemar_preprocess.m | 28 ++++++++++++++++------------
 eeg_nemar_vis.m        | 12 +++++++++++-
 2 files changed, 27 insertions(+), 13 deletions(-)

diff --git a/eeg_nemar_preprocess.m b/eeg_nemar_preprocess.m
index 34859dc..48dd973 100644
--- a/eeg_nemar_preprocess.m
+++ b/eeg_nemar_preprocess.m
@@ -59,6 +59,9 @@
     disp(status_tbl)
     status = table2array(status_tbl);
 
+    splitted = split(EEG.filename(1:end-4),'_');
+    modality = splitted{end};
+
     fprintf('Running pipeline sequence %s\n', strjoin(opt.pipeline, '->'));
     try
         for i=1:numel(opt.pipeline)
@@ -88,18 +91,18 @@
                 rm_chan_types = {'AUDIO','MEG','EOG','ECG','EMG','EYEGAZE','GSR','HEOG','MISC','PPG','PUPIL','REF','RESP','SYSCLOCK','TEMP','TRIG','VEOG'};
                 if isfield(EEG.chanlocs, 'type')
                     EEG = pop_select(EEG, 'rmchantype', rm_chan_types);
-                    types = {EEG.chanlocs.type};
-                    eeg_indices = strmatch('EEG', types)';
-                    if ~isempty(eeg_indices)
-                        EEG = pop_select(EEG, 'chantype', 'EEG');
-                    else
-                        warning("No EEG channel type detected (for first EEG file). Keeping all channels");
-                    end
+		    if strcmp(modality, 'eeg')
+                        types = {EEG.chanlocs.type};
+                        eeg_indices = strmatch('EEG', types)';
+                        if ~isempty(eeg_indices)
+                            EEG = pop_select(EEG, 'chantype', 'EEG');
+                        else
+                            warning("No EEG channel type detected (for first EEG file). Keeping all channels");
+                        end
+	    	    end
                 else
-                    warning("Channel type not detected (for first EEG file)");
+                    warning("Channel type not detected (for first recording file)");
                 end
-                % ALLEEG = pop_select( ALLEEG,'nochannel',{'VEOG', 'Misc', 'ECG', 'M2'});
-
                 status_tbl.remove_chan = 1;
             end
 
@@ -120,12 +123,12 @@
                     'LineNoiseCriterion',4,'Highpass',[0.75 1.25] ,'BurstCriterion',20, ...
                     'WindowCriterion',0.25,'BurstRejection','on','Distance','Euclidian', ...
                     'WindowCriterionTolerances',[-Inf 7] ,'fusechanrej',1}; % based on Arnaud paper
-                % ALLEEG = parexec(ALLEEG, 'pop_clean_rawdata', opt.logdir, options{:});
                 EEG = pop_clean_rawdata( EEG, options{:});
 
                 status_tbl.cleanraw = 1;
             end
 
+	    %{
             if strcmp(operation, "avg_ref")
                 if resume && status_tbl.avg_ref
                     fprintf('Skipping avg_ref\n');
@@ -138,6 +141,7 @@
 
                 status_tbl.avg_ref = 1;
             end
+	    %}
 
             if strcmp(operation, "runica")
                 if resume && status_tbl.runica
@@ -154,7 +158,7 @@
                 status_tbl.runica = 1;
             end
 
-            if strcmp(operation, "iclabel")
+            if strcmp(operation, "iclabel") && strcmp(modality, 'eeg')
                 if resume && status_tbl.iclabel
                     fprintf('Skipping iclabel\n');
                     continue
diff --git a/eeg_nemar_vis.m b/eeg_nemar_vis.m
index cd6a78e..a7a5a22 100644
--- a/eeg_nemar_vis.m
+++ b/eeg_nemar_vis.m
@@ -56,6 +56,9 @@
 
     fprintf('Plots: %s\n', strjoin(opt.plots, ', '));
 
+    splitted = split(EEG.filename(1:end-4), '_');
+    modality = splitted{end};
+
     try
         for i=1:numel(opt.plots)
             plot = opt.plots{i};
@@ -72,7 +75,7 @@
                 plot_IC_activation(EEG);
             end
 
-            if strcmp(plot, 'icmap')
+            if strcmp(plot, 'icmap') && ~strcmp(modality, 'ieeg')
                 plot_ICLabel(EEG);
             end
 
@@ -151,6 +154,10 @@ function plot_spectra(EEG, varargin)
         g = finputcheck(varargin, { 'freq'    'integer' []         [6, 10, 22]; ...
                         'freqrange'   'integer'   []         [1 70]; ...
                         'percent'   'integer'    [], 10});
+	
+	% average reference before plotting
+        EEG = pop_reref(EEG,[], 'interpchan', []);
+
         % spectopo plot
         [spec, freqs] = spectopo(EEG.data, 0, EEG.srate, 'freqrange', g.freqrange, 'title', '', 'chanlocs', EEG.chanlocs, 'percent', g.percent,'plot', 'off');
         [~,ind50]=min(abs(freqs-50));
@@ -180,6 +187,9 @@ function plot_IC_activation(EEG)
             error('No IC decomposition found for EEG')
         end
 
+	% average reference before plotting
+        EEG = pop_reref(EEG,[], 'interpchan', []);
+
         EEG = pop_icflag(EEG,[0.75 1;NaN NaN;NaN NaN;NaN NaN;NaN NaN;NaN NaN;NaN NaN]);
         % IC activations plot
         iclocs = EEG.chanlocs;