make_cortical

function status=make_cortical(sigma_v,sigma_m,bad_chan,bad_obs,infname,outfname)

For now just estimates cortical surface dipoles which is the same thing as a 
Laplacian.  Will eventually integrate evidence etc. once i know how
to use it

This function accepts either 2, 4 or 6 input arguments. If the
last two arguments are omitted the function is run in 
interactive mode, and file names are obtained through the std
dialog boxes. 

sigma_v and sigma_m are variance (noise) estimates on the data and the 
dipole moments respectively

Uses: ave_hdr_offsets_v, rd_egis_hdr_v, zeromeans, rd_onetr_allch

0001 function status=make_cortical(sigma_v,sigma_m,bad_chan,bad_obs,infname,outfname)
0002 %function status=make_cortical(sigma_v,sigma_m,bad_chan,bad_obs,infname,outfname)
0003 %
0004 %For now just estimates cortical surface dipoles which is the same thing as a
0005 %Laplacian.  Will eventually integrate evidence etc. once i know how
0006 %to use it
0007 %
0008 %This function accepts either 2, 4 or 6 input arguments. If the
0009 %last two arguments are omitted the function is run in
0010 %interactive mode, and file names are obtained through the std
0011 %dialog boxes.
0012 %
0013 %sigma_v and sigma_m are variance (noise) estimates on the data and the
0014 %dipole moments respectively
0015 %
0016 %Uses: ave_hdr_offsets_v, rd_egis_hdr_v, zeromeans, rd_onetr_allch
0017 
0018 %Modification History
0019 %    9/95 Created by Jerry and Ramesh
0020 %    3/22/96 bad_obs filtering by BCR
0021 
0022 status = -1;
0023 
0024 %Check number of input arguments
0025 if ~((nargin == 6)|(nargin == 2)|(nargin == 4))
0026     error('Number of input arguments must be 2, 4 or 6.');
0027 end
0028 if nargin == 2
0029     bad_chan = [];bad_obs=[];
0030 end;
0031 if sigma_m == 0
0032     sigma_m = 1;
0033 end;
0034 %Initialize fids
0035 srcfid = -1;
0036 destfid = -1;
0037 %First try batch mode
0038 if nargin == 6
0039     srcfid = fopen(infname, 'r');
0040     if srcfid == -1
0041         error(['Could not open input file' infname '.']);
0042     end
0043     destfid = fopen(outfname, 'wb');
0044     if destfid == -1
0045         temp =fclose(srcfid);
0046         error(['Could not open output file ' outfname '.']);
0047     end
0048 %Otherwise run interactive mode
0049 else
0050     while srcfid == -1
0051         [srcfid, infname]=get_fid('r','*.ave*', 'Open Average File:');
0052     end
0053     while destfid == -1
0054         [destfid, outfname]=put_fid('wb','new.ave_crt','Save New CorticalAverage File As:');
0055     end
0056 end
0057 
0058 %Call EGIS hdr index script
0059 ave_hdr_offsets_v;
0060 %Read in the EGIS header
0061 [fhdr,chdr,ename,czeros,cgains,cnames,fcom,ftext, coff]=rd_egis_hdr_v(srcfid);
0062 %Calculate nyquist frequency and compare to lowpass parameter
0063 % Copy source header to destination file
0064 frewind(srcfid);
0065 [temp, countin] = fread(srcfid, fhdr(LHeader), 'char');
0066 countout = fwrite(destfid, temp, 'char');
0067 if countin ~= countout
0068   error('Error copying header information');
0069 end
0070 if (nargin > 2)&(bad_chan ~= [])
0071     if size(bad_chan,1) ~= chdr(1,NObs)
0072         error('mismatch between bad_chan and NObs')
0073     end;
0074 end;
0075 [xelec,yelec,zelec] = electrodes(fhdr(NChan));
0076 sources = xyz2tp(xelec,yelec,zelec);
0077 A = transfer_matrix(500,0.2,80,8.0,8.2,8.7,9.2,7.8,sources,sources);
0078 
0079 %Begin looping through cells
0080 for c=1:fhdr(NCells)
0081     disp(['Processing cell: ' int2str(c)]);
0082     %Loop through subject averages
0083     for t=1:chdr(c,NObs)
0084         if ~any(t==bad_obs)
0085             trialdata = rd_onetr_allch(srcfid, coff(c), t, fhdr(NChan), chdr(c, NPoints));
0086             trial_mask = ones(1,fhdr(NChan));
0087             if (nargin > 2)&(bad_chan ~= [])
0088                 trial_mask(1,bad_chan(t,find(bad_chan(t,:)))) = zeros(1,size(find(bad_chan(t,:)),2));
0089             end;
0090             good_chan = find(trial_mask);
0091             good_trial = zeros(size(trialdata,1),size(good_chan,2));
0092             good_trial = trialdata(:,good_chan);
0093             good_A = zeros(size(good_chan,2),size(good_chan,2));
0094             good_A = A(good_chan,good_chan);
0095             
0096             
0097             cortical_trial = zeros(size(trialdata,1),size(trialdata,2));
0098             [m,deviations] = bayes_dipole_trial(good_A,good_trial',sigma_v,sigma_m);
0099             cortical_trial(:,good_chan) = m';
0100             cortical_trial = cortical_trial*10;
0101             fwrite(destfid,cortical_trial', 'int16');
0102         else
0103             cortical_trial=zeros(chdr(c, NPoints),fhdr(NChan));
0104             fwrite(destfid,cortical_trial', 'int16');
0105         end
0106     end     %for t=1:chdr(c,NObs)
0107 end     %for c=1:fhdr(NCells)
0108 
0109 disp('Welcome to the cortical surface.  Transfer function fee of $5 applies');
0110 temp=fclose(srcfid);
0111 temp=fclose(destfid);
0112 status = 1;