filtering and down sampling
[data] = vb_filter_raw_data(data,parm)
--- Input
data : input signal [Nchannel x Tsample] or [Ntrial x Tsample]
--- preprocess parameter
parm.freq : Sampling frequency [Hz]
parm
.bias_flg : Bias correction flag (= 0/1/2 : OFF/Bias/Linear)
= N > 5: Bias correction by first N sample
.common_flg : Common reference flag (= 1/0 : ON/OFF)
.highpass : Highpass filter cutoff frequency [Hz]
.lowpass : Lowpass filter cutoff frequency [Hz]
.notch : notch filter stop band frequency = [flow fhigh] [Hz]
.fsamp : Down sampling frequency [Hz]
--- If these fields are empty,
corresponding process are not applied
--- Type of filter specification
.highpass_online : filter order of IIR highpass filter
= 0: FIR highpass filter (eegfilt)
= 1: online highpass filter (exponential)
> 1: Butterworth highpass filter (filtfilt)
< 0: Butterworth highpass filter (online)
.lowpass_online : filter order of IIR lowpass filter
= 0: FIR lowpass filter (eegfilt)
= 1: online lowpass filter (exponential)
> 1: Butterworth lowpass filter (filtfilt)
< 0: Butterworth lowpass filter (online)
2008-5-26 Masa-aki Sato
2009-2-16 Masa-aki Sato
2009-8-10 Masa-aki Sato
Copyright (C) 2011, ATR All Rights Reserved.
License : New BSD License(see VBMEG_LICENSE.txt)0001 function [data] = vb_filter_raw_data(data,parm) 0002 % filtering and down sampling 0003 % [data] = vb_filter_raw_data(data,parm) 0004 % --- Input 0005 % data : input signal [Nchannel x Tsample] or [Ntrial x Tsample] 0006 % --- preprocess parameter 0007 % parm.freq : Sampling frequency [Hz] 0008 % parm 0009 % .bias_flg : Bias correction flag (= 0/1/2 : OFF/Bias/Linear) 0010 % = N > 5: Bias correction by first N sample 0011 % .common_flg : Common reference flag (= 1/0 : ON/OFF) 0012 % .highpass : Highpass filter cutoff frequency [Hz] 0013 % .lowpass : Lowpass filter cutoff frequency [Hz] 0014 % .notch : notch filter stop band frequency = [flow fhigh] [Hz] 0015 % .fsamp : Down sampling frequency [Hz] 0016 % --- If these fields are empty, 0017 % corresponding process are not applied 0018 % --- Type of filter specification 0019 % .highpass_online : filter order of IIR highpass filter 0020 % = 0: FIR highpass filter (eegfilt) 0021 % = 1: online highpass filter (exponential) 0022 % > 1: Butterworth highpass filter (filtfilt) 0023 % < 0: Butterworth highpass filter (online) 0024 % 0025 % .lowpass_online : filter order of IIR lowpass filter 0026 % = 0: FIR lowpass filter (eegfilt) 0027 % = 1: online lowpass filter (exponential) 0028 % > 1: Butterworth lowpass filter (filtfilt) 0029 % < 0: Butterworth lowpass filter (online) 0030 % 0031 % 2008-5-26 Masa-aki Sato 0032 % 2009-2-16 Masa-aki Sato 0033 % 2009-8-10 Masa-aki Sato 0034 % 0035 % Copyright (C) 2011, ATR All Rights Reserved. 0036 % License : New BSD License(see VBMEG_LICENSE.txt) 0037 0038 if (size(data,1) > 1) && (size(data,2)==1), 0039 data = data'; 0040 end 0041 0042 [Nch ,T] = size(data); 0043 0044 fsamp = parm.freq; 0045 0046 % Bias correction 0047 if isfield(parm,'bias_flg') && ~isempty(parm.bias_flg) 0048 if parm.bias_flg > 5, 0049 % Bias correction by first N sample 0050 %fprintf('Bias correction by first %d sample\n',parm.bias_flg) 0051 bias = mean(data(:,1:parm.bias_flg) , 2); 0052 data = vb_repadd(data, - bias); 0053 elseif parm.bias_flg > 0 0054 % Bias correction by all data 0055 %fprintf('Bias correction mode = %d\n',parm.bias_flg) 0056 data = vb_channel_bias_remove(data,parm.bias_flg); 0057 end 0058 end 0059 0060 % Online highpass filter 0061 if isfield(parm,'highpass') && ~isempty(parm.highpass), 0062 if isfield(parm,'highpass_online') 0063 % Highpass cutoff frequency 0064 fcut = parm.highpass; 0065 Norder = parm.highpass_online; 0066 0067 if Norder == 1 0068 % Exponential online highpass filter 0069 %fprintf('Online highpass filter %6.1f [Hz]\n',fcut) 0070 data = vb_online_highpass_cut(data, fsamp, fcut); 0071 elseif Norder > 1 0072 %fprintf('butter filtfilt highpass filter %6.1f [Hz]\n',fcut) 0073 [B,A] = butter(Norder, fcut/(fsamp/2) ,'high'); 0074 for n=1:Nch 0075 data(n,:) = filtfilt(B,A, data(n,:) ); 0076 end 0077 elseif Norder < 0 0078 %fprintf('Online butter highpass filter %6.1f [Hz]\n',fcut) 0079 [A, B, C, D, Z] = ... 0080 vb_online_butter_init(fcut,fsamp,abs(Norder),Nch,'high'); 0081 if Nch > 10 0082 [data ,Z] = vb_online_filter_loop(A, B, C, D, data, Z); 0083 else 0084 [data ,Z] = vb_online_filter_loop_mex(A, B, C, D, data, Z); 0085 end 0086 end 0087 end 0088 end 0089 0090 % Lowpass filter 0091 if isfield(parm,'lowpass') && ~isempty(parm.lowpass) 0092 fcut = parm.lowpass; 0093 0094 if isfield(parm,'lowpass_online') 0095 % fprintf('Online Lowpass filter %6.1f [Hz]\n', fcut) 0096 % Online Butterworth filter 0097 Norder = parm.lowpass_online; % filter order 0098 if Norder > 0 0099 %fprintf('butter filtfilt lowpass filter %6.1f [Hz]\n',fcut) 0100 [B,A] = butter(Norder, fcut/(fsamp/2) ,'low'); 0101 for n=1:Nch 0102 data(n,:) = filtfilt(B,A, data(n,:) ); 0103 end 0104 elseif Norder < 0 0105 %fprintf('Online butter lowpass filter %6.1f [Hz]\n',fcut) 0106 [A, B, C, D, Z] = ... 0107 vb_online_butter_init(fcut,fsamp,abs(Norder),Nch,'low'); 0108 if Nch > 10 0109 [data ,Z] = vb_online_filter_loop(A, B, C, D, data, Z); 0110 else 0111 [data ,Z] = vb_online_filter_loop_mex(A, B, C, D, data, Z); 0112 end 0113 end 0114 end 0115 if ~isfield(parm,'lowpass_online') || parm.lowpass_online==0 0116 % FIR Lowpass filter 0117 % fprintf('FIR Lowpass filter %6.1f [Hz]\n', fcut) 0118 data = eegfilt(data, fsamp, 0, fcut ); 0119 end 0120 end 0121 0122 % IIR notch filter 0123 if isfield(parm,'notch') && ~isempty(parm.notch), 0124 Norder = 5; 0125 [B,A] = butter(Norder, fcut/(fsamp/2) ,'stop'); 0126 0127 for m=1:Ntr 0128 for n=1:Nch 0129 data(n,:,m) = filtfilt(B,A, data(n,:,m) ); 0130 end 0131 end 0132 end 0133 0134 % Down sampling 0135 if isfield(parm,'fsamp') && ~isempty(parm.fsamp) && fsamp ~=parm.fsamp 0136 fsnew = parm.fsamp; 0137 %fprintf('Down sampling to %6.1f [Hz]\n',fsnew) 0138 0139 data = vb_convert_freq(data, fsamp, fsnew); 0140 fsamp = fsnew; 0141 end 0142 0143 % FIR Highpass filter 0144 if isfield(parm,'highpass') && ~isempty(parm.highpass) 0145 if ~isfield(parm,'highpass_online') || parm.highpass_online==0 0146 %fprintf('FIR Highpass filter %6.1f [Hz]\n',parm.highpass) 0147 data = eegfilt(data, fsamp, parm.highpass, 0); 0148 end 0149 end 0150 0151 % Common reference 0152 if isfield(parm,'common_flg') && ~isempty(parm.common_flg) && parm.common_flg==1 0153 % fprintf('Common reference\n') 0154 data = vb_repadd(data, - mean(data,1)); 0155 end