set2avg() - Computes ERPs from an EEGLAB set file or EEG variable and writes them to a Kutaslab average file. ERPs will be normalized via the cal pulse ERP and can be artifact corrected via ICA (if artifact ICs have been labeled). Thus the average file produced is analogous to a Kutaslab ".nrm" file. Usage: >> set2avg(EEG_or_setname,out_fname,hdr_fname,varargin) Required Inputs: EEG_or_setname - [EEGLAB EEG variable | EEGLAB set filename] Identifies the EEGLAB EEG variable in memory or the set file on disk from which the ERPs will be generated. If specifying a filename, be sure to include the file's path unless the file is in the current working directory. IC's corresponding to EEG artifacts should be labelled as such in the field EEG.iclabels. See the function is_art.m for a list of IC labels recognized as artifacts. out_fname - [string] The root of the name of the file(s) that will be created (e.g., 'sbjct1'). The file extension ".nrm", ".icnrm", ".iccav" will be automatically added. ".nrm" is used for ERPs derived from only artifact free trials (no ICA artifact correction). ".icnrm" is the same as ".nrm" but with ICA artifact correction (useful for seeing how ICA might be distorting ERPs). ".iccor" is used for ERPs derived from all trials with ICA artifact correction. hdr_fname - [string] The name of an existing Kutaslab average file (e.g., sbjct2.nrm) or a header file produced by the program mdh. The header from this file will be used as a template for the new file that is being created. Any average or header file should work (i.e., it does not need to have the same channels or bins as the set file). The template file cannot have the same name as one of the average files being produced by set2avg. Optional Inputs: 'bsln' - [start_time stop_time | 'allpre'] Time window (in milliseconds) that will be used to baseline each trial. In other words, for a given trial and electrode the mean voltage in that time window will be removed from each time point in the trial. [] means no baselining (data will be averaged as it is in the set file after ICA artifact correction). 'allpre' means that all pre-event time points (including Time 0) will be used for baselining. Note, if you perform ICA artifact correction, you'll want to baseline the data here as any previous baselining might no longer hold. {Default: []} 'pp10uv' - [scalar] Data resolution ("points per microvolt") at which the averages will be produced. This parameter determines how the data are scaled before writing to the average file. Ideally, you want to make the data as big as possible, while remaining in the limits of numeric representation range. Kutaslab default for this is 1000. {Default: 1000} 'ica_correct' - [1 | 0] If 1, ICA will be used to correct for EEG artifacts and the average file produced will have a ".iccor" extension. ICs are identified as artifacts via the field EEG.iclabels. See is_art.m for a list of labels recognized as artifacts. If 0, ICA is not used to correct for artifacts and the average file produced will have a ".nrm" extension. {Default: 0} 'just_cor_drty' - [1 | 0] If 1, ICA artifact correction will only be applied to polluted trials and ERPs will be formed from both artifact-free and corrected trials. This argument has no effect if the argument 'ica_correct' is set to 0. {Default: 0} 'cln_avgs' - [1 | 0] If 1, two additional average files will be produced with the same root filename as out_fname. The first of these will have the extension ".nrm" and will be formed only from artifact free trials (artifact free according to the arf file used to import the data into EEGLAB). The second of these will have the extension ".icnrm" and will be composed of the exact same ERPs in the ".nrm" file after the ICA artifact correction filter has been applied. Contrasting these two average files is useful for estimating how ICA artifact correction might be distorting the data. {Default: 0} 'overwrite' - [1 | 0] If 1, existing average files that have the same name as the new files being produced will be automatically overwritten. If 0, the user will be asked if s/he wants to overwrite them. {Default: 0} 'verblevel' - an integer specifiying the amount of information you want functions to provide about what they are doing during runtime. Options are: 0 - quiet, only show errors, warnings, and EEGLAB reports 1 - stuff anyone should probably know 2 - stuff you should know the first time you start working with a data set {default value if not already globally specified} 3 - stuff that might help you debug (show all reports) query_ccode - [1 | 0] If one, the user will be asked at the command line for the condition code of bins with zero trials (if any). If zero, the function will simply use the condition code for the previous bin. {default: 0} Outputs: An average file with a name specified by 'out_fname' is written to disk. Nothing is output to MATLAB. If 'cln_avgs' option is used, two additional average files will be produced. Example: >> set_fname='/homes/dgroppe/SANDBOX/PRSENT/prsent41.set'; >> hdr_fname='/homes/dgroppe/SANDBOX/PRSENT/prsent.hdr'; >> set2avg(set_fname,'prsent41',hdr_fname,'bsln',[-100 0],'cln_avgs',1,'ica_correct',1); Author: David Groppe Kutaslab, 11/2009 Notes: -This function relies on erpio to create the avg files. Thus it only works on the Kutaslab network (unless you have a compiled copy of erpio on whatever machine you're using). -Code was written to minimize memory demands at the cost of being slower. If faster code is needed, the function could be revised.
0001 function set2avg(EEG_or_setname,out_fname,hdr_fname,varargin) 0002 % set2avg() - Computes ERPs from an EEGLAB set file or EEG variable and writes 0003 % them to a Kutaslab average file. ERPs will be normalized 0004 % via the cal pulse ERP and can be artifact corrected via ICA 0005 % (if artifact ICs have been labeled). Thus the average file 0006 % produced is analogous to a Kutaslab ".nrm" file. 0007 % 0008 % Usage: 0009 % >> set2avg(EEG_or_setname,out_fname,hdr_fname,varargin) 0010 % 0011 % Required Inputs: 0012 % EEG_or_setname - [EEGLAB EEG variable | EEGLAB set filename] Identifies 0013 % the EEGLAB EEG variable in memory or the set file on disk 0014 % from which the ERPs will be generated. If specifying a 0015 % filename, be sure to include the file's path unless the 0016 % file is in the current working directory. IC's 0017 % corresponding to EEG artifacts should be labelled as 0018 % such in the field EEG.iclabels. See the function 0019 % is_art.m for a list of IC labels recognized as artifacts. 0020 % out_fname - [string] The root of the name of the file(s) that will be 0021 % created (e.g., 'sbjct1'). The file extension ".nrm", 0022 % ".icnrm", ".iccav" will be automatically added. ".nrm" 0023 % is used for ERPs derived from only artifact free trials 0024 % (no ICA artifact correction). ".icnrm" is the same as 0025 % ".nrm" but with ICA artifact correction (useful for 0026 % seeing how ICA might be distorting ERPs). ".iccor" is 0027 % used for ERPs derived from all trials with ICA artifact 0028 % correction. 0029 % hdr_fname - [string] The name of an existing Kutaslab average file 0030 % (e.g., sbjct2.nrm) or a header file produced by the program 0031 % mdh. The header from this file will be used as a 0032 % template for the new file that is being created. Any 0033 % average or header file should work (i.e., it does not 0034 % need to have the same channels or bins as the set 0035 % file). The template file cannot have the same name as 0036 % one of the average files being produced by set2avg. 0037 % 0038 % Optional Inputs: 0039 % 'bsln' - [start_time stop_time | 'allpre'] Time window (in milliseconds) 0040 % that will be used to baseline each trial. In other words, for 0041 % a given trial and electrode the mean voltage in that time window 0042 % will be removed from each time point in the trial. [] means no 0043 % baselining (data will be averaged as it is in the set file after ICA 0044 % artifact correction). 'allpre' means that all pre-event time 0045 % points (including Time 0) will be used for baselining. Note, if 0046 % you perform ICA artifact correction, you'll want to baseline the 0047 % data here as any previous baselining might no longer hold. 0048 % {Default: []} 0049 % 'pp10uv' - [scalar] Data resolution ("points per microvolt") at which the 0050 % averages will be produced. This parameter determines how the 0051 % data are scaled before writing to the average file. Ideally, 0052 % you want to make the data as big as possible, while remaining 0053 % in the limits of numeric representation range. Kutaslab default 0054 % for this is 1000. {Default: 1000} 0055 % 'ica_correct' - [1 | 0] If 1, ICA will be used to correct for EEG artifacts 0056 % and the average file produced will have a ".iccor" 0057 % extension. ICs are identified as artifacts via the field 0058 % EEG.iclabels. See is_art.m for a list of labels recognized 0059 % as artifacts. If 0, ICA is not used to correct for artifacts 0060 % and the average file produced will have a ".nrm" extension. 0061 % {Default: 0} 0062 % 'just_cor_drty' - [1 | 0] If 1, ICA artifact correction will only be 0063 % applied to polluted trials and ERPs will be formed from both 0064 % artifact-free and corrected trials. This argument has no 0065 % effect if the argument 'ica_correct' is set to 0. {Default: 0} 0066 % 'cln_avgs' - [1 | 0] If 1, two additional average files will be produced 0067 % with the same root filename as out_fname. The first of 0068 % these will have the extension ".nrm" and will be formed 0069 % only from artifact free trials (artifact free according 0070 % to the arf file used to import the data into EEGLAB). 0071 % The second of these will have the extension ".icnrm" 0072 % and will be composed of the exact same ERPs in the ".nrm" file 0073 % after the ICA artifact correction filter has been applied. 0074 % Contrasting these two average files is useful for estimating 0075 % how ICA artifact correction might be distorting the data. 0076 % {Default: 0} 0077 % 'overwrite' - [1 | 0] If 1, existing average files that have the same 0078 % name as the new files being produced will be automatically 0079 % overwritten. If 0, the user will be asked if s/he wants 0080 % to overwrite them. {Default: 0} 0081 % 'verblevel' - an integer specifiying the amount of information you want 0082 % functions to provide about what they are doing during runtime. 0083 % Options are: 0084 % 0 - quiet, only show errors, warnings, and EEGLAB reports 0085 % 1 - stuff anyone should probably know 0086 % 2 - stuff you should know the first time you start working 0087 % with a data set {default value if not already globally 0088 % specified} 0089 % 3 - stuff that might help you debug (show all reports) 0090 % query_ccode - [1 | 0] If one, the user will be asked at the 0091 % command line for the condition code of bins with zero 0092 % trials (if any). If zero, the function will simply 0093 % use the condition code for the previous bin. {default: 0} 0094 % 0095 % 0096 % Outputs: 0097 % An average file with a name specified by 'out_fname' is written to disk. 0098 % Nothing is output to MATLAB. If 'cln_avgs' option is used, two 0099 % additional average files will be produced. 0100 % 0101 % 0102 % Example: 0103 % >> set_fname='/homes/dgroppe/SANDBOX/PRSENT/prsent41.set'; 0104 % >> hdr_fname='/homes/dgroppe/SANDBOX/PRSENT/prsent.hdr'; 0105 % >> set2avg(set_fname,'prsent41',hdr_fname,'bsln',[-100 0],'cln_avgs',1,'ica_correct',1); 0106 % 0107 % Author: 0108 % David Groppe 0109 % Kutaslab, 11/2009 0110 % 0111 % Notes: 0112 % -This function relies on erpio to create the avg files. Thus it only 0113 % works on the Kutaslab network (unless you have a compiled copy of erpio 0114 % on whatever machine you're using). 0115 % 0116 % -Code was written to minimize memory demands at the cost of being slower. 0117 % If faster code is needed, the function could be revised. 0118 0119 0120 %%%%%%%%%%%%%%%% REVISION LOG %%%%%%%%%%%%%%%%% 0121 % 0122 %12/11/09-Added ic_distortion option to measure degree of ERP distortion 0123 %caused by ICA artifact correction. Also, code now uses inputParser. -DG 0124 % 0125 %2/10/09-Made code able to deal with EEG.epoch(trial).eventccode being a cell 0126 %array (which it is if you run eeg_checkset.m on it). -DG 0127 % 0128 %7/29/10-Added query_ccode option and made it possible for set2avg to deal 0129 %with bins that have no trials. Also added a check to make sure that the 0130 %number of channels of cal pulse information matches the number of channels 0131 %in the data.-DG 0132 0133 %Input Parser 0134 p=inputParser; 0135 %Note: the order of the required arguments needs to match precisely their 0136 %order in the function definition (which is also the order used by p.parse 0137 %below) 0138 p.addRequired('EEG_or_setname',@(x) isstruct(x) || ischar(x)); 0139 p.addRequired('out_fname',@ischar); 0140 p.addRequired('hdr_fname',@ischar); 0141 p.addParamValue('est_ic_dstrt',0,@isnumeric); 0142 p.addParamValue('dstrt_bins',[],@isnumeric); %default is to use all bins 0143 p.addParamValue('bsln','allpre',@(x) isnumeric(x) || strcmpi(x,'allpre')); %default is to use all prestim time points 0144 p.addParamValue('pp10uv',1000,@isnumeric); %default is to use all time points 0145 p.addParamValue('verblevel',[],@(x) x>=0); 0146 p.addParamValue('cln_avgs',0,@isnumeric); 0147 p.addParamValue('just_cor_drty',0,@isnumeric); 0148 p.addParamValue('ica_correct',0,@isnumeric); 0149 p.addParamValue('overwrite',0,@isnumeric); 0150 p.addParamValue('query_ccode',0,@isnumeric); 0151 0152 p.parse(EEG_or_setname,out_fname,hdr_fname,varargin{:}); 0153 0154 global VERBLEVEL 0155 0156 if isempty(VERBLEVEL) && isempty(p.Results.verblevel) 0157 VERBLEVEL=2; 0158 elseif ~isempty(p.Results.verblevel) 0159 VERBLEVEL=p.Results.verblevel; 0160 end 0161 0162 % Load/manage EEG struct 0163 if ischar(EEG_or_setname), 0164 %EEG_or_setname is a set filename 0165 EEG=pop_loadset(EEG_or_setname); %note this reports an informative error if the set file doesn't exist 0166 % It also reports on the command line what file it's attempting to load 0167 else 0168 %EEG_or_setname is an EEGLAB EEG struct 0169 EEG=EEG_or_setname; %note that this does not create a second copy of EEG 0170 %in memory since neither variable is yet modified. 0171 clear EEG_or_setname; 0172 end 0173 0174 n_bin=length(EEG.bindesc); 0175 [n_chan, n_tpt, n_epoch]=size(EEG.data); 0176 0177 %Check to make sure number of cal pulse channels equals the number of 0178 %EEG.data channels 0179 n_cal_pulse=size(EEG.cal_info.erps,1); 0180 if n_cal_pulse~=n_chan, 0181 error('You have %d channels of cal pulse averages (in EEG.cal_info.erps) but %d channels of EEG data (in EEG.data).\n', ... 0182 n_cal_pulse,n_chan); 0183 end 0184 n_cal_pulse=size(EEG.cal_info.usedcal_itemnum,1); 0185 if n_cal_pulse~=n_chan, 0186 error('You have %d channels of cal pulse item numbers (in EEG.cal_info.usedcal_itemnum) but %d channels of EEG data (in EEG.data).\n', ... 0187 n_cal_pulse,n_chan); 0188 end 0189 0190 % Check to make sure header template exists 0191 [status, w]=unix(['ls ' hdr_fname]); 0192 if status, 0193 error('Header template file %s doesn''t exist.',hdr_fname); 0194 end 0195 0196 % Check to make sure number of time points is a multiple of 256 0197 rmndr=rem(n_tpt,256); 0198 if rmndr, 0199 %fprintf('The number of time points per trial is not a multiple of 256.\n'); 0200 n_pad=256-rmndr; 0201 padding=zeros(n_chan,n_pad); 0202 %fprintf('Adding %d time points of zeros to the end of each trial to correct this.\n',n_pad); 0203 VerbReport(sprintf(['The number of time points per trial is not a multiple of 256.\n' ... 0204 'Adding %d time points of zeros to the end of each trial to correct this.\n'],n_pad),1,VERBLEVEL); 0205 else 0206 n_pad=0; 0207 padding=[]; 0208 end 0209 0210 0211 % Get baseline window 0212 if ~isempty(p.Results.bsln), 0213 if strcmpi(p.Results.bsln,'allpre'), 0214 %use all prestimulus time points 0215 bsln_tpt(1)=1; 0216 bsln_tpt(2)=find_tpt(0,EEG.times); 0217 else 0218 bsln_tpt(1)=find_tpt(p.Results.bsln(1),EEG.times); 0219 bsln_tpt(2)=find_tpt(p.Results.bsln(2),EEG.times); 0220 end 0221 bsln_win=bsln_tpt(1):bsln_tpt(2); 0222 if isempty(bsln_win), 0223 error('Baseline window from %.1f to %.1f ms includes no time points. Please revise it.', ... 0224 bsln_tpt(1),bsln_tpt(2)); 0225 else 0226 VerbReport(sprintf('Removing mean EEG from %.1f ms to %.1f ms for each epoch/trial.\n', ... 0227 EEG.times(bsln_tpt(1)),EEG.times(bsln_tpt(2))),1,VERBLEVEL); 0228 end 0229 end 0230 0231 0232 0233 0234 %% Remove artifact ICs if they've been labeled 0235 0236 % Check to see if some independent components have been labeled as artifacts 0237 fldnames=fieldnames(EEG); 0238 ic_arts=0; 0239 for a=1:length(fldnames), 0240 if strcmpi(fldnames{a},'iclabels'), 0241 ic_arts=1; 0242 end 0243 end 0244 0245 if ic_arts, 0246 n=length(EEG.iclabels); 0247 art_ics=[]; 0248 for dg=1:n, 0249 if is_art(EEG.iclabels{dg}), 0250 art_ics=[art_ics dg]; 0251 end 0252 end 0253 0254 unmix=EEG.icaweights*EEG.icasphere; 0255 mix=EEG.icawinv; 0256 0257 if p.Results.est_ic_dstrt, 0258 crnt_level=VERBLEVEL; 0259 %Set VERBLEVEL 1 to suppress redudant report about baselining data 0260 ic_distortion(EEG,'bins',p.Results.dstrt_bins,'bsln',p.Results.bsln,'plot',1,'verblevel',1); 0261 VERBLEVEL=crnt_level; % Reset VERBLEVEL to what it was before 0262 drawnow; 0263 else 0264 if VERBLEVEL>=1, 0265 fprintf('%d ICs have been labeled as EEG artifacts.\n',length(art_ics)); 0266 if ~isempty(art_ics), 0267 fprintf('The ICs labeled by artifacts are: \n'); 0268 for dg=art_ics, 0269 fprintf('IC %d: %s\n',dg,EEG.iclabels{dg}); 0270 end 0271 end 0272 end 0273 eeglab_rej=find(EEG.reject.gcompreject==1); 0274 dif1=setdiff(eeglab_rej,art_ics); 0275 dif2=setdiff(art_ics,eeglab_rej); 0276 if ~isempty(dif1) || ~isempty(dif2) 0277 if isempty(eeglab_rej), 0278 rej_ic_str='None'; 0279 else 0280 rej_ic_str=int2str(eeglab_rej); 0281 end 0282 msg=['The ICs marked for rejection by EEGLAB in EEG.reject.gcompreject differ from those labeled as artifacts in EEG.iclabels.' ...); 0283 10 'Only the artifact labels in EEG.iclabels will be used.' 10 ... 0284 'The ICs labeled by artifacts by EEG.reject.gcompreject are: ' rej_ic_str]; 0285 watchit(msg); 0286 end 0287 end 0288 0289 if p.Results.cln_avgs, 0290 % Create an avg file from just clean trials (i.e., clean according to arf 0291 % file) 0292 just_clean_trials=1; 0293 ica_correct=0; 0294 VerbReport('Making ERPs from just artifact free trials, without ICA artifact correction.',2,VERBLEVEL); 0295 write_avg(EEG,out_fname,hdr_fname,bsln_win,p.Results.pp10uv,just_clean_trials,ica_correct,p.Results.overwrite,n_pad,padding,art_ics,0,p.Results.query_ccode); 0296 0297 % Filter the ERPs from just clean trials with ICA and write another average file. 0298 % Useful for seeing how much ICA might be distorting the data 0299 just_clean_trials=1; 0300 ica_correct=1; 0301 VerbReport('Making ERPs from just artifact free trials, with ICA artifact correction.',2,VERBLEVEL); 0302 write_avg(EEG,out_fname,hdr_fname,bsln_win,p.Results.pp10uv,just_clean_trials,ica_correct,p.Results.overwrite,n_pad,padding,art_ics,0,p.Results.query_ccode); 0303 end 0304 0305 if p.Results.ica_correct, 0306 if p.Results.just_cor_drty 0307 VerbReport('Making ERPs from all trials, with ICA artifact correction applied only to polluted trials.',2,VERBLEVEL); 0308 else 0309 VerbReport('Making ERPs from all trials, with ICA artifact correction applied to all trials.',2,VERBLEVEL); 0310 end 0311 just_clean_trials=0; 0312 write_avg(EEG,out_fname,hdr_fname,bsln_win,p.Results.pp10uv,just_clean_trials,p.Results.ica_correct,p.Results.overwrite,n_pad,padding,art_ics,p.Results.just_cor_drty,p.Results.query_ccode); 0313 end 0314 else 0315 if p.Results.ica_correct, 0316 error(['This EEG variable does not have ICA components. You need to run ICA first before you can perform ICA ' ... 0317 'artifact correction.']); 0318 else 0319 just_clean_trials=1; 0320 VerbReport('Making ERPs from just clean trials, without ICA artifact correction.',2,VERBLEVEL); 0321 art_ics=[]; 0322 end 0323 write_avg(EEG,out_fname,hdr_fname,bsln_win,p.Results.pp10uv,just_clean_trials,p.Results.ica_correct,p.Results.overwrite,n_pad,padding,art_ics,0,p.Results.query_ccode); 0324 end 0325 0326 0327 0328 function write_avg(EEG,out_fname,hdr_fname,bsln_win,pp10uv,just_clean_trials,ica_correct,ovrwrite,n_pad,padding,art_ics,just_cor_drty,query_ccode) 0329 %function write_avg(EEG,out_fname,hdr_fname,bsln_win,pp10uv,just_clean_trials,ica_correct,ovrwrite,n_pad,padding,art_ics,just_cor_drty,query_ccode) 0330 % 0331 % Inputs: 0332 % EEG - EEGLAB EEG variable from which ERPs will be derived. 0333 % out_fname - Name of average file that will be produce. The following 0334 % naming conventions are automatically tagged on: 0335 % *.drty=all trials, no ICA artifact correction 0336 % *.nrm=just clean trials (according to garv log flags), no 0337 % ICA artifact correction 0338 % *.iccor=all trials, yes ICA artifact correction 0339 % *.icnrm=just clean trials (according to garv log flags), 0340 % yes ICA artifact correction (useful for estimating 0341 % ERP distortion produced by ICA correction) 0342 % hdr_fname - Name of the header file or average file that will be used as 0343 % a template for the new average file being produced. This file 0344 % will not be modified. 0345 % bsln_win - The time window (in milliseconds) that will be used to 0346 % baseline each trial 0347 % pp10uv - Points per 10 microvolts (determines scale of data in avg 0348 % file) 0349 % just_clean_trials - [1 | 0] If 1, only clean trials (according to garv 0350 % log flags) will be used to form ERPs. 0351 % ica_correct - [1 | 0] If 1, ICA will be used to remove EEG 0352 % artifacts from the EEG. 0353 % ovrwrite - [1 | 0] If 1, existing average files that have the 0354 % same name as the new files being produced will be 0355 % automatically overwritten. If 0, user will be 0356 % asked if s/he wants to overwrite them. 0357 % n_pad - Number of 0 microvolt time points that will be 0358 % added to ERPs so that the number of time points are 0359 % a multiple of 256. 0360 % art_ics - [binary vector] If the ith element of art_ics is a 0361 % 1, the ith IC corresponds to an EEG artifact. 0362 % just_cor_drty - [1 | 0] If 1, ICA artifact correction correction 0363 % is only applied to polluted trials (according to garv 0364 % log flags) 0365 % query_ccode - [scalar] If non-zero, the user will be asked at the 0366 % command line for the condition code of bins with zero 0367 % trials (if any). If zero, the function will simply 0368 % use the condition code for the previous bin. 0369 0370 global VERBLEVEL; 0371 0372 % Data parameters 0373 s=size(EEG.data); 0374 n_tpt=s(2); 0375 n_chan=s(1); 0376 n_epoch=s(3); 0377 n_bin=length(EEG.bindesc); 0378 0379 % Create avg file name 0380 if ~just_clean_trials && ica_correct, 0381 ext='.iccor'; 0382 elseif just_clean_trials && ica_correct, 0383 ext='.icnrm'; 0384 elseif ~just_clean_trials && ~ica_correct, 0385 ext='.drty'; % I don't know why would you ever want to make such an average, 0386 % but the option is included just in case. It is not commented though in set2avg.m. 0387 else 0388 % just clean trials but no ica correction 0389 ext='.nrm'; 0390 end 0391 0392 if ~ismember('.',out_fname), 0393 %tack on auto-extension 0394 out_fname=[out_fname ext]; 0395 else 0396 id=find(out_fname=='.'); 0397 mx_id=max(id); 0398 if mx_id==1, 0399 error('Output filename %s cannot begin with a period.',out_fname); 0400 end 0401 %remove existing extension & tack on auto-extension 0402 out_fname=[out_fname(1:(mx_id-1)) ext]; 0403 end 0404 0405 0406 % Check to make sure out_fname doesn't already exist 0407 [status, w]=unix(['ls ' out_fname]); 0408 if ~status, 0409 while ~ovrwrite 0410 rsp=input(sprintf('File %s already exists. Do you want to overwrite it (Y or N)?',out_fname),'s'); 0411 if strcmpi(rsp,'Y') || strcmpi(rsp,'yes'), 0412 ovrwrite=1; 0413 elseif strcmpi(rsp,'N') || strcmpi(rsp,'no'), 0414 ovrwrite=-1; 0415 end 0416 end 0417 0418 if ovrwrite==1, 0419 % Remove current version of output file if it already exists 0420 VerbReport(sprintf('Removing old version of file %s.',out_fname),1,VERBLEVEL); 0421 [status, w]=unix(['rm -f ' out_fname]); % ?? add an error check? 0422 %Note, status and w provide no information about the success of rm -f 0423 elseif ovrwrite==-1, 0424 fprintf('\n'); % line break to make command line reports look better 0425 return; 0426 end 0427 end 0428 0429 0430 % Create avg file 0431 VerbReport(sprintf('Creating average file %s using header information from %s.', ... 0432 out_fname,hdr_fname),1,VERBLEVEL); 0433 fid=erpio('createavg',out_fname,hdr_fname); 0434 if (fid==-1) 0435 error('Could not create ERP file %s using erpio.',out_fname); 0436 end 0437 erpio('put_hdrvar',fid,'pp10uv',pp10uv); 0438 0439 0440 erpio('put_hdrvar',fid,'cprecis',(n_tpt+n_pad)/256); 0441 erpio('put_hdrvar',fid,'expdesc',EEG.comments); 0442 erpio('put_hdrvar',fid,'subdesc',EEG.subject); 0443 erpio('put_hdrvar',fid,'presampling',-EEG.times(1)); 0444 %erpio('put_hdrvar',fid,'epoch',EEG.times(end)-EEG.times(1)+1000/EEG.srate); 0445 %--at least one of my set files had EEG.times in units of seconds (weird) 0446 erpio('put_hdrvar',fid,'epoch',1000*(n_pad+n_tpt)/EEG.srate); 0447 erpio('put_hdrvar',fid,'chans',n_chan); 0448 %erpio('put_hdrvar',fid,'prfdesc','average'); ?? makes matlab crash 0449 erpio('put_hdrvar',fid,'clktick',100000/EEG.srate); 0450 for c=1:n_chan, 0451 erpio('put_hdrvar',fid,'chndesc',c-1,EEG.chanlocs(c).labels); 0452 end 0453 0454 %pre-allocate memory 0455 erps=zeros(n_chan,n_tpt,n_bin); 0456 erp_ct=zeros(1,n_bin); 0457 ccode=-ones(1,n_bin); 0458 0459 % create ICA filter (if requested) 0460 if ica_correct, 0461 non_art=setdiff(1:n_chan,art_ics); 0462 unmix=EEG.icaweights*EEG.icasphere; 0463 mix=inv(unmix); %just in case EEG.icawinv has been scaled differently than mixing matrix 0464 ic_filt=mix(:,non_art)*unmix(non_art,:); 0465 end 0466 0467 %add up all trials 0468 for trial=1:n_epoch, 0469 if iscell(EEG.epoch(trial).eventlogflag) 0470 eventLogFlag=0; 0471 for lfLoop=1:length(EEG.epoch(trial).eventlogflag) 0472 eventLogFlag=eventLogFlag+EEG.epoch(trial).eventlogflag{lfLoop}; 0473 end 0474 else 0475 eventLogFlag=EEG.epoch(trial).eventlogflag; 0476 end 0477 if ~just_clean_trials || ~eventLogFlag 0478 % clean trials have an eventlogflag of 0 0479 if ica_correct, 0480 if just_cor_drty && eventLogFlag, 0481 sweep=ic_filt*squeeze(EEG.data(:,:,trial)); 0482 elseif ~just_cor_drty, 0483 sweep=ic_filt*squeeze(EEG.data(:,:,trial)); 0484 else 0485 sweep=squeeze(EEG.data(:,:,trial)); 0486 end 0487 else 0488 sweep=squeeze(EEG.data(:,:,trial)); 0489 end 0490 sweep=rmbase(sweep,n_tpt,bsln_win); 0491 for types=1:length(EEG.epoch(trial).eventtype), 0492 if strcmpi(EEG.epoch(trial).eventtype{types}(1:3),'bin'), 0493 bin_num=str2num(EEG.epoch(trial).eventtype{types}(4:end)); 0494 erps(:,:,bin_num)=erps(:,:,bin_num)+sweep; 0495 erp_ct(bin_num)=erp_ct(bin_num)+1; 0496 if iscell(EEG.epoch(trial).eventccode), % if eeg_checkset has been run eventccode will be a cell array 0497 if ccode(bin_num)==-1, 0498 ccode(bin_num)=EEG.epoch(trial).eventccode{types}; 0499 elseif ccode(bin_num)~=EEG.epoch(trial).eventccode{types}, 0500 error('Trial %d has a Condition Code of %d and contributes to Bin %d, but at least one other trial in that bin has a Condition Code of %d.', ... 0501 trial,EEG.epoch(trial).eventccode{types},bin_num,ccode(bin_num)); 0502 end 0503 else 0504 if ccode(bin_num)==-1, 0505 ccode(bin_num)=EEG.epoch(trial).eventccode; 0506 elseif ccode(bin_num)~=EEG.epoch(trial).eventccode, 0507 error('Trial %d has a Condition Code of %d and contributes to Bin %d, but at least one other trial in that bin has a Condition Code of %d.', ... 0508 trial,EEG.epoch(trial).eventccode,bin_num,ccode(bin_num)); 0509 end 0510 end 0511 end 0512 end 0513 end 0514 end 0515 0516 crnt_cond=0; 0517 for b=0:n_bin, 0518 if b==0, 0519 erpio('put_hdrvar',fid,'condcode',0); 0520 erpio('put_hdrvar',fid,'condesc','Amp Calibration'); 0521 erpio('put_hdrvar',fid,'sums',EEG.cal_info.npulse_used); 0522 erpio('put_hdrvar',fid,'totrawrecs',EEG.cal_info.npulse_psbl); 0523 erpio('writebin',fid,[EEG.cal_info.erps padding]*pp10uv/10,'Calibration Pulses'); 0524 %note erpio returns nothing for 'writebin' (i.e. no indication 0525 %of success) 0526 else 0527 if ccode(b)<0, 0528 if query_ccode, 0529 fprintf('There are no trials for Bin %d thus I don''t know what its condition code is.\n',b); 0530 ccode(b)=input(['Enter the condition code for Bin ' int2str(b) ': ']); 0531 else 0532 if VERBLEVEL>=2, 0533 fprintf('There are no trials for Bin %d. I will use the condition code for the previous bin (i.e., %d).\n',b,crnt_cond); 0534 end 0535 ccode(b)=crnt_cond; 0536 end 0537 end 0538 if crnt_cond~=ccode(b), 0539 crnt_cond=ccode(b); 0540 erpio('put_hdrvar',fid,'condcode',crnt_cond); 0541 erpio('put_hdrvar',fid,'condesc',EEG.condesc{crnt_cond}); 0542 end 0543 if VERBLEVEL>=2, 0544 fprintf('Bin %d, CCode %d, Trials %d(%d): %s\n',b,ccode(b),erp_ct(b),EEG.rawtrials_per_bin(b),EEG.bindesc{b}); 0545 end 0546 %write avg 0547 erpio('put_hdrvar',fid,'sums',erp_ct(b)); 0548 erpio('put_hdrvar',fid,'totrawrecs',EEG.rawtrials_per_bin(b)); 0549 if erp_ct(b), 0550 erpio('writebin',fid,[squeeze(erps(:,:,b)*pp10uv/(10*erp_ct(b))) padding],EEG.bindesc{b}); 0551 else 0552 erpio('writebin',fid,[squeeze(erps(:,:,b)) padding],EEG.bindesc{b}); 0553 end 0554 end 0555 end 0556 0557 %close file 0558 erpio('close',fid); 0559 0560 fprintf('\n'); % line break to make command line reports look better 0561 0562