ctf_read_res4_coef.m

含有多种ICA算法的eeglab工具箱

if VERBOSE,
  fprintf('...number of filters : %d\n',ctf.res4.numberFilters);
  for i = 1:ctf.res4.numberFilters,
    fprintf(' -Filter # \t%g\n',i)
    fprintf(' -Frequency: \t%g Hz\n',ctf.res4.filter(i).frequency)
    fprintf(' -Class: \t%g\n',ctf.res4.filter(i).class)
    fprintf(' -Type: \t%g\n',ctf.res4.filter(i).type)
    if ~isempty(ctf.res4.filter(i).parameters)
      fprintf(' -Parameter(s): \t%g\n',ctf.res4.filter(i).parameters)
    end
  end
end





%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% READ SENSOR INFORMATION

% seek channel information byte offset
fseek(res4fid,b+f,-1);

% read channel names (not trivial!)
for channel = 1:ctf.res4.generalSetup.numberChannels,
  temp = fread(res4fid,[1,32],'char');
  temp(temp>127) = 0;
  temp(temp<0) = 0;
  temp = strtok(temp,char(0));
  temp = strtok(temp,'-'); % this is in brainstorm code
  ctf.sensor.info(channel).label = char(temp);
  ctf.sensor.label{1,channel} = char(temp);
end


for chan = 1:ctf.res4.generalSetup.numberChannels,
  
  %ftell(res4fid);
  
  ctf.sensor.info(chan).index = fread(res4fid,1,'int16');
  ctf.sensor.info(chan).extra = fread(res4fid,1,'int16');
  id = fread(res4fid,1,'int32')+1;
  ctf.sensor.info(chan).proper_gain = fread(res4fid,1,'double');
  ctf.sensor.info(chan).q_gain = fread(res4fid,1,'double');
  ctf.sensor.info(chan).io_gain = fread(res4fid,1,'double');
  ctf.sensor.info(chan).io_offset = fread(res4fid,1,'double');
  fread(res4fid,1,'int16');
  fread(res4fid,1,'int16');
  fread(res4fid,1,'int32');
  
  %fseek(res4fid,ftell(res4fid)+6,0);
  
  for pos = 1:8,
    ctf.sensor.info(chan).coil(pos).position.x = fread(res4fid,1,'double');
    ctf.sensor.info(chan).coil(pos).position.y = fread(res4fid,1,'double');
    ctf.sensor.info(chan).coil(pos).position.z = fread(res4fid,1,'double');
    fread(res4fid,1,'double');
    ctf.sensor.info(chan).coil(pos).orient.x = fread(res4fid,1,'double');
    ctf.sensor.info(chan).coil(pos).orient.y = fread(res4fid,1,'double');
    ctf.sensor.info(chan).coil(pos).orient.z = fread(res4fid,1,'double');
    fread(res4fid,1,'double');
    fread(res4fid,1,'int16');
    fread(res4fid,1,'int32');
    fread(res4fid,1,'int16');
    fread(res4fid,1,'double');
    
    %fseek(res4fid,ftell(res4fid)+56,0);
    %fseek(res4fid,ftell(res4fid)-80,0);
  end
  
  for pos = 1:8,
    ctf.sensor.info(chan).hcoil(pos).position.x = fread(res4fid,1,'double');
    ctf.sensor.info(chan).hcoil(pos).position.y = fread(res4fid,1,'double');
    ctf.sensor.info(chan).hcoil(pos).position.z = fread(res4fid,1,'double');
    fread(res4fid,1,'double');
    ctf.sensor.info(chan).hcoil(pos).orient.x = fread(res4fid,1,'double');
    ctf.sensor.info(chan).hcoil(pos).orient.y = fread(res4fid,1,'double');
    ctf.sensor.info(chan).hcoil(pos).orient.z = fread(res4fid,1,'double');
    fread(res4fid,1,'double');
    fread(res4fid,1,'int16');
    fread(res4fid,1,'int32');
    fread(res4fid,1,'int16');
    fread(res4fid,1,'double');
    
    %fseek(res4fid,ftell(res4fid)+56,0);
    %fseek(res4fid,ftell(res4fid)+80,0);
  end
  %fseek(res4fid,ftell(res4fid)+1288,-1);
end




%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Find channel types

ctf.sensor.type.eeg = 9;
ctf.sensor.type.meg = 5;
ctf.sensor.type.ref = [0 1];


%ctf.sensor.type.vc  = ??;
%vcsens = find([ctf.sensor.info.index] == ??);


eegsens = find([ctf.sensor.info.index] == 9);
megsens = find([ctf.sensor.info.index] == 5);
refsens = find([ctf.sensor.info.index] == 1);
refsens = [refsens,find([ctf.sensor.info.index] == 0)];
allsens = [eegsens,megsens,refsens];
othersens = setdiff([1:ctf.res4.generalSetup.numberChannels],allsens);

ctf.sensor.index = struct(...
  'eeg',{eegsens},...
  'meg',{megsens},...
  'ref',{refsens},...
  'other',{othersens} ...
  );



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Location coordinates of channels in centimeters, in patient head space

for chan = 1:ctf.res4.generalSetup.numberChannels,
  
  switch ctf.sensor.info(chan).index,
    
    case {0,1,5},
      %0=Reference Channels, 
      %1=More Reference Channels, 
      %5=MEG Channels
      
      coord = [ctf.sensor.info(chan).hcoil(1:2).position];
      ctf.sensor.info(chan).location = [coord.x; coord.y; coord.z];
      
      orient = [ctf.sensor.info(chan).hcoil(1).orient];
      ctf.sensor.info(chan).orientation = [orient.x; orient.y; orient.z];
      
      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      % This ensures that the orientation of the sensor is away from the
      % center of a sphere.  It uses the sign of the dot product between
      % the orientation vector and the location vector.
      tmp = ctf.sensor.info(chan).orientation' * ctf.sensor.info(chan).location;
      tmp = sign(tmp(1));
      ctf.sensor.info(chan).orientation = tmp * ctf.sensor.info(chan).orientation;
      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      
    case 9,
      %EEG Channels    
      
      coord = [ctf.sensor.info(chan).hcoil(1:2).position];
      ctf.sensor.info(chan).location = [coord.x; coord.y; coord.z];
      ctf.sensor.info(chan).orientation = [];
      
  end
end








%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% NUMBER OF COEFFICIENTS, byte offset = b+f+nc*1360, byte size = 1

offset = b + f + ctf.res4.generalSetup.numberChannels * 1360;
fseek(res4fid,offset,-1);
ctf.res4.numberCoefficients = fread(res4fid,1,'int16');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% SENSOR COEFFICIENT RECORDS, byte offset = b+f+nc*1360+2, byte size = 1992

if VERBOSE & ctf.res4.numberCoefficients,
  fprintf('...reading %d coefficients\n',ctf.res4.numberCoefficients);
end

SENSOR_LABEL  = 31;
MAX_NUM_COEFS = 50;
MAX_BALANCING = MAX_NUM_COEFS;

hexadef = {'00000000','47314252','47324252','47334252','47324f49','47334f49'};
strdef = {'NOGRAD','G1BR','G2BR','G3BR','G2OI','G3OI'};

for i = 1:ctf.res4.numberCoefficients,
  
  % read the sensor name (channel label)
  temp = fread(res4fid,[1,32],'char');
  temp(temp>127) = 0;
  temp(temp<0) = 0;
  temp = strtok(temp,char(0));
  temp = strtok(temp,'-');
  sensorName = char(temp);
  sensorIndex = strmatch( sensorName, ctf.sensor.label );
  
  % read the coefficient type
  coefType = fread(res4fid,1,'bit32');
  
  padding = fread(res4fid,1,'int32'); % not sure why this is needed???
  
  % read the coefficient record
  numberCoefs = fread(res4fid,1,'int16');
  
  if numberCoefs > MAX_NUM_COEFS,
    msg = sprintf('numberCoefs > MAX_NUM_COEFS\n');
    warning(msg);
  end
  
  sensor_list = char(fread(res4fid,[SENSOR_LABEL,MAX_BALANCING],'uchar')');
  % clean-up the sensor_list
  sensor_list = sensor_list(1:numberCoefs,:);
  for i=1:numberCoefs,
    
    temp = strtok(sensor_list(i,:),char(0));
    temp = strtok(temp,'-');
    % check if this sensor is a reference
    refLabels = ctf.sensor.label(refsens);
    refIndex = strmatch(temp,refLabels);
    if refIndex,
      % ensure this one has the same label
      temp = refLabels{refIndex};
    end
    new_sensor_list(i,:) = temp;
    
  end
  sensor_list = new_sensor_list;
  
  coefs_list = fread(res4fid,MAX_BALANCING,'double');
  % clean-up the coefs_list
  coefs_list = coefs_list(1:numberCoefs,:);
  
  % allocate the coefficient parameters into the ctf struct
  ctf.res4.sensorCoef(i).sensorName = sensorName;
  ctf.res4.sensorCoef(i).coefType   = coefType;
  ctf.res4.sensorCoef(i).coefRec.numberCoefs = numberCoefs;
  ctf.res4.sensorCoef(i).coefRec.sensor_list = sensor_list;
  ctf.res4.sensorCoef(i).coefRec.coefs_list  = coefs_list;
  
  
  
%   % DLW:
%   % This is a brainstorm variable, note the use of coefType
%   % Not clear why this is checked and allocated as such
%   % It is not returned by this function, as yet
%   coefType = find( hex2dec(hexadef) == coefType );
%   if coefType,
%     CoefInfo{sensorIndex,coefType-1}.numberCoefs = numberCoefs;
%     for i=1:numberCoefs,
%       CoefInfo{sensorIndex,coefType-1}.sensor_list(i) = sensor_list(i);
%       CoefInfo{sensorIndex,coefType-1}.coefs(i) = coefs_list(i);
%     end
%   end
  
  
  
end


t = toc; fprintf('...done (%6.2f sec)\n\n',t);

return