interp_by_samp

function status = interp_by_samp(cells,obs,samps,plot_res,interp_type,bad_chan,output_root,header_flag,infname, stack)

cells = rows by columns array of cell mapping to page
obs = rows by columns array of obs mapping to page
samps = vector of samples to interpolate
plot_res = grid size
interp_type - 'sp' = spherical harmonic = '3d' for 3-d splines
 bad_chan = array or matrix of bad_channels 
 output_root = root name for output image
 header_Flag = flag top determine swhether the images get a header, 'h' for header 's' for no header defaults of 'h'
 stack = if this is true, interpolations will be placed into 1 file

0001 function status = interp_by_samp(cells,obs,samps,plot_res,interp_type,bad_chan,output_root,header_flag,infname,stack)
0002 %function status = interp_by_samp(cells,obs,samps,plot_res,interp_type,bad_chan,output_root,header_flag,infname, stack)
0003 %
0004 %cells = rows by columns array of cell mapping to page
0005 %obs = rows by columns array of obs mapping to page
0006 %samps = vector of samples to interpolate
0007 %plot_res = grid size
0008 %interp_type - 'sp' = spherical harmonic = '3d' for 3-d splines
0009 % bad_chan = array or matrix of bad_channels
0010 % output_root = root name for output image
0011 % header_Flag = flag top determine swhether the images get a header, 'h' for header 's' for no header defaults of 'h'
0012 % stack = if this is true, interpolations will be placed into 1 file
0013 
0014 %Modification History
0015 %    2/95 Created by Ramesh Srinivasan
0016 %          included features from PJ's update to make_interp_by_samp which is now defunct
0017 %  02/23/96  PJ  nobs variable is set equal to number of rows in obs.  There is
0018 %                still a fundamental problem with the way obs and bad_chan
0019 %                are handled.
0020 %
0021 %  01/21/97 PJ Added option to put all interpolations into a single file
0022 %
0023 %  04/26/97 PJ Updated to handle more efficient routines in /spline
0024 %              Handling of more than one obs is almost certainly defunct.
0025 %
0026 %  11/23/97 PJ Fixed masking problem under Matlab 5.0
0027 
0028 %Check number of input arguments
0029 if ~((nargin ==9)|(nargin == 5)|(nargin == 6)|(nargin == 8)|(nargin == 10))
0030     error('Number of input arguments must be either 5 6 8 or 9 or 10.');
0031 end
0032 %Initialize fids
0033 srcfid = -1;
0034 
0035 if ~((size(obs,1) == size(cells,1))&(size(obs,2) == size(cells,2)))
0036     error('page layout is faulty')
0037 end;
0038 if nargin == 5
0039     bad_chan = [];
0040     header_flag = 'h';
0041     output_root = 'interpolation';
0042 end;
0043 if isempty(output_root)
0044     output_root = [infname '_interp_'];
0045 end;
0046 if ~((header_flag == 'h')|(header_flag == 's'))
0047     header_flag = 'h';
0048 end;
0049 if isempty(samps),  samps = 1:chdr(1:NPoints), end;
0050 
0051 if nargin < 10
0052   stack = 0;
0053 end
0054 
0055 %First try batch mode
0056 if nargin == 9 | nargin == 10
0057     srcfid = fopen(infname, 'r','ieee-be');
0058     if srcfid == -1
0059         error(['Could not open input file' infname '.']);
0060     end
0061 
0062 %Otherwise run interactive mode
0063 else
0064     while srcfid == -1
0065         [srcfid, infname]=get_fid('r','*.*', 'Open File:');
0066     end
0067 end
0068 
0069 %Call EGIS hdr index script
0070 ave_hdr_offsets_v;
0071 
0072 %Read in the EGIS header
0073 [fhdr,chdr,ename,czeros,cgains,cnames,fcom,ftext, coff]=rd_egis_hdr_v(srcfid);
0074 
0075 % CREATE CHANNEL MASK
0076 good_chan_mask = ones(chdr(1,NObs),fhdr(NChan));
0077 
0078 nobs = size(obs,1);
0079 
0080 if ~isempty(bad_chan)
0081   % Check to see if number of rows in obs and bad_chans jives
0082   nrows_bad = size(bad_chan,1);
0083   if nobs ~= nrows_bad
0084     disp(['WARNING: Number of observations (' int2str(nobs) ') doesn''t match number of rows in bad_chan (' int2str(nrows_bad) ')']);
0085     if nrows_bad == 1
0086       disp('Replicating bad_chan vector to match nobs');
0087       bad_chan = ones(size(good_chan_mask,1))*bad_chan;
0088     else
0089       error('Drastic mismatch between nobs and bad_chan')
0090     end
0091   end
0092 
0093   for isub = 1:size(good_chan_mask,1)
0094     bchan = find(bad_chan(isub,:));
0095     good_chan_mask(isub,bad_chan(isub,bchan)) = zeros(1,size(bchan,2));
0096   end;
0097 end;
0098 
0099 [xelec,yelec,zelec] = electrodes(fhdr(NChan));
0100 if interp_type == 'sp'
0101     if fhdr(NChan) == 129
0102         nmax = 8;
0103     else
0104         error('something is rotten')
0105     end;
0106 end;
0107 
0108 good_chan = ones(1,fhdr(NChan));
0109 if ~isempty(bad_chan)
0110   good_chan(bad_chan) = zeros(size(bad_chan));
0111 end
0112 
0113 good_chan = find(good_chan);
0114 
0115 [xelec,yelec,zelec] = electrodes(129);
0116 
0117   v = zeros(1,size(good_chan,2));
0118   x = zeros(1,size(good_chan,2));
0119   y = zeros(1,size(good_chan,2));
0120   z = zeros(1,size(good_chan,2));
0121 
0122   x = xelec(good_chan);
0123   y = yelec(good_chan);
0124   z = zelec(good_chan);
0125 
0126 welec = 1;
0127 [k, kinv, a, ainv, e] = k_and_e(welec,x,y,z);
0128 
0129 % WORK THRHOUGH EACH SAMPLE
0130 
0131 for isamp = 1:size(samps,2)
0132   disp(['Working on sample: ' int2str(samps(isamp))]);
0133   image_binary = zeros(plot_res*size(cells,1),plot_res*size(cells,2));
0134   for icell=1:size(cells,1)
0135         for t=1:size(cells,2)
0136         c = cells(icell,t);
0137 %        disp(['  Cell: ' int2str(c)]);
0138 %        disp(['    Obs: ' int2str(obs(icell,t))]);
0139            trialdata = rd_onetr_allch(srcfid, coff(c), obs(icell,t), fhdr(NChan), chdr(c, NPoints));
0140         w = trialdata(samps(isamp),:);
0141         v = w(good_chan);
0142 
0143         if interp_type == 'sp'
0144             meanv = mean(v);
0145             v = v - meanv*ones(1,size(v,2));
0146               [sp_coff,error_check, sp_mat] = bayes_sphcoff(nmax,x,y,z,v);
0147               imagemat = image_2D(sp_coff,nmax,plot_res,120,meanv);
0148         elseif interp_type == '3d'
0149               [p,q,error_check]= mateqs(welec,x,y,z,v,k,kinv,a,ainv,e);
0150             ruu = sqrt(x(1).^2+y(1).^2+z(1).^2);
0151             [xs,ys,zs] = polgrid(plot_res,ruu);    
0152             image_3d= interp_3d(welec,x,y,z,xs,ys,zs,p,q);    
0153             imagemat= mask(xs,ys,zs,image_3d,120); 
0154         else
0155             error('havent dealt with spherical splines yet')
0156         end;            
0157         image_binary(plot_res*(icell-1)+1:plot_res*(icell-1)+plot_res,plot_res*(t-1)+1:plot_res*(t-1)+plot_res)=imagemat;
0158          
0159     end;  %for t=1:size(obs,2)
0160   end     %for c=1:size(cells,2)
0161 
0162   if ~stack  
0163     outfname = [output_root '_s' int2str(samps(isamp))];
0164     outfid = fopen(outfname,'wb','ieee-be');
0165   elseif isamp == 1
0166     outfname = [output_root '_s' int2str(samps(isamp)) '_s' int2str(samps(length(samps)))]
0167     outfid = fopen(outfname,'wb','ieee-be');
0168   end
0169   
0170   if isamp == 1 | ~stack
0171     if header_flag == 'h'    
0172       num_written = fwrite(outfid, size(image_binary), 'float');
0173         if num_written ~= 2
0174         error('Failed to write header info');
0175         end
0176 
0177       fwrite(outfid, plot_res, 'float');
0178     end;
0179   end
0180       num_written = fwrite(outfid,image_binary','float');
0181       if num_written ~= size(image_binary,1)*size(image_binary,2)
0182         error('write-out filed')
0183       end;
0184   if ~stack
0185     fclose(outfid);
0186   end
0187 end     %for c=1:size(samps,2)
0188 
0189 disp('Finished running make_interp.');
0190 temp=fclose('all');
0191 status = 1;