eidors_readdata

PURPOSE ^

EIDORS readdata - read data files from various EIT equipment

SYNOPSIS ^

function [vv, auxdata, stim ]= eidors_readdata( fname, format, frame_range, extra )

DESCRIPTION ^

 EIDORS readdata - read data files from various EIT equipment
    manufacturers

 [vv, auxdata, stim ]= eidors_readdata( fname, format, frame_range, extra )
    frame_range is the list of frames to load from the file (ie. 1:100).
    if the frames are beyond the number in the file, no data is returned
    to get all frames, use frame_range= []

 Currently the list of supported file formats is:
    - "EIT" Files, may be one of 
          - Draeger Pulmovista (2008+)
          - GoeIIMF/Carefusion (2010+)
          - Swisstom BB2/Pioneer (2010+)
       The function will attempt to autodetect the format
       
    - MCEIT (GoeIIMF / Viasys) "get" file format 
        format = "GET" or "MCEIT"
        Note that the data is "untwisted" to correspond to a "no_rotate_meas" stim pattern
    - MCEIT (GoeIIMF) "get" file format
        format = "GET-RAW"
        Data in original order, corresponds to "rotate_meas" stim pattern

    - Draeger (Pulmovista) "EIT" file format (2008+)
        format = "DRAEGER-EIT"

    - Draeger "get" file format (- 2007 - format for Draeger equipment)
        format = "DRAEGER-GET"

    - Swisstom EIT equipment "EIT" (Pioneer set and BB2):
           'LQ1' (2010 - 2011)
           'LQ2' (2013 - 2014)
           'LQ4' (2015+)

    - Dixtal file format, from Dixtal inc, Brazil
        format = 'DIXTAL_encode' extract encoder from provided Dll
           - Output: [encodepage] = eidors_readdata( path,'DIXTAL_encode');
              where path= '/path/to/Criptografa_New.dll' (provided with system)
        format = 'DIXTAL'
           - output: [vv] = eidors_readdata( fname, 'DIXTAL', [], encodepage );
    - New Carefusion "EIT" file format
        format = "GOEIIMF-EIT" or "carefusion"

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function [vv, auxdata, stim ]= eidors_readdata( fname, format, frame_range, extra )
0002 % EIDORS readdata - read data files from various EIT equipment
0003 %    manufacturers
0004 %
0005 % [vv, auxdata, stim ]= eidors_readdata( fname, format, frame_range, extra )
0006 %    frame_range is the list of frames to load from the file (ie. 1:100).
0007 %    if the frames are beyond the number in the file, no data is returned
0008 %    to get all frames, use frame_range= []
0009 %
0010 % Currently the list of supported file formats is:
0011 %    - "EIT" Files, may be one of
0012 %          - Draeger Pulmovista (2008+)
0013 %          - GoeIIMF/Carefusion (2010+)
0014 %          - Swisstom BB2/Pioneer (2010+)
0015 %       The function will attempt to autodetect the format
0016 %
0017 %    - MCEIT (GoeIIMF / Viasys) "get" file format
0018 %        format = "GET" or "MCEIT"
0019 %        Note that the data is "untwisted" to correspond to a "no_rotate_meas" stim pattern
0020 %    - MCEIT (GoeIIMF) "get" file format
0021 %        format = "GET-RAW"
0022 %        Data in original order, corresponds to "rotate_meas" stim pattern
0023 %
0024 %    - Draeger (Pulmovista) "EIT" file format (2008+)
0025 %        format = "DRAEGER-EIT"
0026 %
0027 %    - Draeger "get" file format (- 2007 - format for Draeger equipment)
0028 %        format = "DRAEGER-GET"
0029 %
0030 %    - Swisstom EIT equipment "EIT" (Pioneer set and BB2):
0031 %           'LQ1' (2010 - 2011)
0032 %           'LQ2' (2013 - 2014)
0033 %           'LQ4' (2015+)
0034 %
0035 %    - Dixtal file format, from Dixtal inc, Brazil
0036 %        format = 'DIXTAL_encode' extract encoder from provided Dll
0037 %           - Output: [encodepage] = eidors_readdata( path,'DIXTAL_encode');
0038 %              where path= '/path/to/Criptografa_New.dll' (provided with system)
0039 %        format = 'DIXTAL'
0040 %           - output: [vv] = eidors_readdata( fname, 'DIXTAL', [], encodepage );
0041 %    - New Carefusion "EIT" file format
0042 %        format = "GOEIIMF-EIT" or "carefusion"
0043 
0044 %    - Sheffield MK I "RAW" file format
0045 %        format = "RAW" or "sheffield"
0046 %    - ITS (International Tomography Systems)
0047 %        format = "ITS" or "p2k"
0048 %    - IIRC (Impedance Imaging Research Center, Korea)
0049 %        format = "txt" or "IIRC"
0050 %    - University of Cape Town formats
0051 %        format = "UCT_SEQ"  UCT sequence file
0052 %           - Output: [ stimulation, meas_select]= eidors_readdata(fname, 'UCT_SEQ')
0053 %        format = "UCT_CAL"  UCT calibration file
0054 %           - Output: [vv, no_cur_caldata_raw ]= eidors_readdata( fname, 'UCT_CAL' )
0055 %                 where no_cur_caldata_raw is data captured with no current
0056 %        format = "UCT_DATA"  UCT data frame file
0057 %           - Output: [vv]= eidors_readdata( fname, 'UCT_DATA' )
0058 %
0059 % Usage
0060 % [vv, auxdata, stim ]= eidors_readdata( fname, format )
0061 %     vv      = measurements - data frames in each column
0062 %     auxdata = auxillary data - if provided by system
0063 %     stim    = stimulation structure, to be used with
0064 %                fwd_model.stimulation.
0065 %     fname = file name
0066 %     stim.framerate = acquisition rate (frames/sec) if available
0067 %
0068 %  if format is unspecified, an attempt to autodetect is made
0069 
0070 % (C) 2005-09 Andy Adler. License: GPL version 2 or version 3
0071 % $Id: eidors_readdata.m 4987 2015-05-12 18:00:21Z aadler $
0072 
0073 % TODO:
0074 %   - output an eidors data object
0075 %   - test whether file format matches specified stimulation pattern
0076 %   - todo MCEIT provides curr and volt on driven electrodes.
0077 %       how can this be provided to system?
0078 
0079 if ~exist(fname,'file')
0080    error([fname,' does not exist']);
0081 end
0082 
0083 if nargin < 2
0084 % unspecified file format, autodetect
0085    dotpos = find(fname == '.');
0086    if isempty( dotpos ) 
0087       error('file format unspecified, can`t autodetect');
0088    else
0089       dotpos= dotpos(end);
0090       format= fname( dotpos+1:end );
0091    end
0092 end
0093 
0094 
0095 auxdata = []; % default, can be overriden if the format has it
0096 fmt = pre_proc_spec_fmt( format, fname );
0097 switch fmt
0098    case 'mceit';
0099       [vv,curr,volt,auxdata_out] = mceit_readdata( fname );
0100       auxdata.auxdata = auxdata_out;
0101       auxdata.curr    = curr;
0102       auxdata.volt    = volt;
0103 
0104       stim = basic_stim(16);
0105    case 'draeger-get'
0106       vv = draeger_get_readdata( fname );
0107 
0108       stim = basic_stim(16);
0109 
0110    case 'draeger-eit'
0111      [fr] = read_draeger_header( fname );
0112      % Currently Draeger equipment uses this pattern with 5mA injection
0113      stim = mk_stim_patterns(16,1,[0,1],[0,1],{'no_rotate_meas','no_meas_current'},.005);
0114      [stim(:).framerate] = deal(fr);
0115      [vv] = read_draeger_file( fname );
0116      auxdata = vv; % the actual voltages
0117 %    vv = vv(1:208,:) + 1j*vv(323+(1:208),:); % I THINK THAT's WHAT THIS IS
0118      vv = vv(1:208,:);
0119 
0120    case {'raw', 'sheffield'}
0121       vv = sheffield_readdata( fname );
0122 
0123       stim = basic_stim(16);
0124    case {'p2k', 'its'}
0125       vv = its_readdata( fname );
0126 
0127       stim = 'UNKNOWN';
0128    case {'txt','iirc'}
0129       vv = iirc_readdata( fname );
0130 
0131       stim = basic_stim(16);
0132    case 'uct_seq'
0133       [vv,auxdata] = UCT_sequence_file( fname );
0134 
0135       stim = 'UNKNOWN';
0136    case 'uct_cal'
0137       [vv,auxdata] = UCT_calibration_file( fname );
0138 
0139       stim = 'UNKNOWN';
0140    case 'uct_data'
0141       [vv] = UCT_LoadDataFrame( fname );
0142 
0143       stim = 'UNKNOWN';
0144    case {'carefusion','goeiimf-eit'}
0145       [vv] = carefusion_eit_readdata( fname );
0146   
0147       stim = basic_stim(16);
0148 
0149    case 'lq1'
0150       [vv] = landquart1_readdata( fname );
0151 
0152       stim = mk_stim_patterns(32,1,[0,5],[0,5],{'no_rotate_meas','no_meas_current'},.005);
0153       
0154    case {'lq2','lq3'}
0155       [vv] = landquart2_readdata( fname );
0156 
0157       stim = mk_stim_patterns(32,1,[0,5],[0,5],{'no_rotate_meas','no_meas_current'},.005);
0158      
0159    case 'lq4'
0160       [vv] = landquart4_readdata( fname );
0161 
0162       stim = mk_stim_patterns(32,1,[0,5],[0,5],{'no_rotate_meas','no_meas_current'},.005);
0163       
0164    case 'dixtal_encode'
0165       [vv] = dixtal_read_codepage( fname );
0166       stim = 'N/A';
0167 
0168    case 'dixtal'
0169       [vv] = dixtal_read_data( fname, frame_range, extra );
0170       auxdata = vv(1025:end,:);
0171       vv      = vv([1:1024],:);
0172       stim= mk_stim_patterns(32,1,[0,4],[0,4], ... 
0173               {'meas_current','no_rotate_meas'}, 1);
0174 
0175    otherwise
0176       error('eidors_readdata: file "%s" format unknown', format);
0177 end
0178 
0179 function stim = basic_stim(N_el);
0180    stim= mk_stim_patterns(16,1,[0,1],[0,1], ... \
0181           {'no_meas_current','no_rotate_meas'}, 1);
0182 
0183 function fmt = pre_proc_spec_fmt( format, fname );
0184    fmt= lower(format);
0185    if strcmp(fmt,'get')
0186       if is_get_file_a_draeger_file( fname)
0187          fmt= 'draeger-get';
0188       else
0189          fmt= 'mceit';
0190       end
0191    end
0192 
0193    if strcmp(fmt,'eit')
0194       draeger =   is_eit_file_a_draeger_file( fname );
0195       swisstom=   is_eit_file_a_swisstom_file( fname );
0196       carefusion= is_eit_file_a_carefusion_file( fname );
0197       if carefusion
0198          eidors_msg('"%s" appears to be in GOEIIMF/Carefusion format',fname,3);
0199          fmt= 'carefusion';
0200       elseif draeger
0201          eidors_msg('"%s" appears to be in Draeger format',fname,3);
0202          fmt= 'draeger-eit';
0203       elseif swisstom
0204          fmt= sprintf('lq%d',swisstom);
0205          eidors_msg('"%s" appears to be in %s format',fname,upper(fmt),3);
0206       else
0207          error('EIT file specified, but it doesn''t seem to be a Carefusion file')
0208       end
0209    end
0210 
0211 function df= is_get_file_a_draeger_file( fname)
0212    fid= fopen(fname,'rb');
0213    d= fread(fid,[1 26],'uchar');
0214    fclose(fid);
0215    df = all(d == '---Draeger EIT-Software---');
0216 
0217 function df= is_eit_file_a_draeger_file( fname );
0218    fid= fopen(fname,'rb');
0219    d= fread(fid,[1 80],'uchar');
0220    fclose(fid);
0221    ff = findstr(d, '---Draeger EIT-Software');
0222    if ff;
0223       df = 1;
0224       eidors_msg('Draeger format: %s', d(ff(1)+(0:30)),4);
0225    else 
0226       df = 0;
0227    end
0228 
0229 function df= is_eit_file_a_swisstom_file( fname );
0230    fid= fopen(fname,'rb');
0231    d= fread(fid,[1 6],'uchar');
0232    fclose(fid);
0233    
0234    if any(d(4)==[2,3,4]) && all(d([1,2,3,5,6]) == 0);
0235       df = d(4);
0236       eidors_msg('Swisstom format: LQ%d', df, 4);
0237    else 
0238       df = 0;
0239    end
0240 
0241 function df = is_eit_file_a_carefusion_file( fname )
0242    fid= fopen(fname,'rb');
0243    d= fread(fid,[1 80],'uchar');
0244    fclose(fid);
0245    df = 0;
0246    if d(1:2) ~= [255,254]; return; end
0247    if d(4:2:end) ~= 0; return; end
0248    str= setstr(d(3:2:end));
0249    tests1= '<?xml version="1.0" ?>';
0250    tests2= '<EIT_File>';
0251    if ~any(findstr(str,tests1)); return; end
0252    if ~any(findstr(str,tests2)); return; end
0253    
0254    df= 1;
0255    return
0256 
0257 function [vv] = carefusion_eit_readdata( fname );
0258    fid= fopen(fname,'rb');
0259    d= fread(fid,[1 180],'uchar');
0260    str= setstr(d(3:2:end));
0261    outv = regexp(str,'<Header>(\d+) bytes</Header>','tokens');
0262    if length(outv) ~=1; 
0263       error('format problem reading carefusion eit files');
0264    end
0265    header = str2num(outv{1}{1});
0266 
0267 % NOTE: we're throwing away the header completely. This isn't
0268 % the 'right' way to read a file.
0269 
0270    fseek(fid, header, -1); % From the beginning
0271    d_EIT= fread(fid,[inf],'float32');
0272 
0273    fseek(fid, header, -1); % From the beginning
0274    d_struct= fread(fid,[inf],'int32');  % In the struct, meaning of every int: 1 Type, 3 sequence No., 4 size in byte, 5 count
0275    fclose(fid);
0276    pt_EIT=1;               % pointer of the EIT data
0277    vv=[];
0278    while (pt_EIT<=length(d_struct))
0279       switch d_struct(pt_EIT)
0280          case 3, % type=3,  EIT transimpedance measurement
0281            % d_struct(pt_EIT+2), Sequence No., start from 0
0282            vv(1:256,1+d_struct(pt_EIT+2))= d_EIT(pt_EIT+6:pt_EIT+6+255);
0283            pt_EIT=pt_EIT+6+256;
0284          case 7, %type=7, EIT image vector
0285            pt_EIT=pt_EIT+912+6;        % drop the image
0286          case 8, %type=8, Waveform data
0287            %drop
0288            pt_EIT=pt_EIT+6+d_struct(pt_EIT+3)/4*d_struct(pt_EIT+4);
0289          case 10,% type = 10, unknown type
0290            %drop
0291            pt_EIT=pt_EIT+6+d_struct(pt_EIT+3)/4*d_struct(pt_EIT+4);
0292          otherwise,
0293            eidors_msg('WARNING: unknown type in carefusion file type');
0294            pt_EIT=pt_EIT+6+d_struct(pt_EIT+3)/4*d_struct(pt_EIT+4);
0295        end
0296    end
0297    vv=vv(47+[1:208],:);
0298    
0299 
0300 % Read the old <2006 Draeger "get" file format
0301 function [vv,curr,volt,auxdata] = draeger_get_readdata( fname );
0302    fid= fopen(fname,'rb');
0303    emptyctr=0;
0304    while emptyctr<2
0305      line= fgetl(fid);
0306      if isempty(line)
0307         emptyctr= emptyctr+1;
0308      else
0309         eidors_msg(line,0);
0310         emptyctr=0;
0311      end
0312    end
0313    d= fread(fid,inf,'float',0,'ieee-le');
0314    fclose(fid);
0315 
0316    if rem( length(d), 256) ~=0
0317       eidors_msg('File length strange - cropping file',0);
0318       d=d(1:floor(length(d)/256)*256);
0319    end
0320 
0321    dd= reshape( d, 256, length(d)/256);
0322    vv= untwist(dd);
0323 
0324    curr=0.00512*dd(209:224,:);  % Amps
0325    volt=12*dd(225:240,:); %Vrms
0326    auxdata= dd(241:255,:);
0327    auxdata= auxdata(:);
0328     
0329 function jnk
0330 %[adler@adler01 sept07]$ xxd Sch_Pneumoperitoneum_01_001.get | head -30
0331 %0000000: 2d2d 2d44 7261 6567 6572 2045 4954 2d53  ---Draeger EIT-S
0332 %0000010: 6f66 7477 6172 652d 2d2d 0d0a 2d2d 2d50  oftware---..---P
0333 %0000020: 726f 746f 636f 6c20 4461 7461 2d2d 2d2d  rotocol Data----
0334 %0000030: 2d0d 0a0d 0a44 6174 653a 2020 3138 2d30  -....Date:  18-0
0335 %0000040: 322d 3230 3034 0d0a 5469 6d65 3a20 2031  2-2004..Time:  1
0336 %0000050: 333a 3138 2050 4d0d 0a0d 0a46 696c 656e  3:18 PM....Filen
0337 %0000060: 616d 653a 2020 2020 2020 2020 2053 6368  ame:         Sch
0338 %0000070: 5f50 6e65 756d 6f70 6572 6974 6f6e 6575  _Pneumoperitoneu
0339 %0000080: 6d5f 3031 5f30 3031 2e67 6574 0d0a 4453  m_01_001.get..DS
0340 %0000090: 5020 5365 7269 616c 204e 722e 3a20 2020  P Serial Nr.:
0341 %00000a0: 4549 5430 322f 3035 2f30 3030 360d 0a0d  EIT02/05/0006...
0342 %00000b0: 0a46 7265 7175 656e 6379 205b 487a 5d3a  .Frequency [Hz]:
0343 %00000c0: 2020 2020 3937 3635 362e 330d 0a47 6169      97656.3..Gai
0344 %00000d0: 6e3a 2020 2020 2020 2020 2020 2020 2020  n:
0345 %00000e0: 2020 2031 3134 350d 0a41 6d70 6c69 7475     1145..Amplitu
0346 %00000f0: 6465 3a20 2020 2020 2020 2020 2020 2031  de:            1
0347 %0000100: 3030 300d 0a53 616d 706c 6520 5261 7465  000..Sample Rate
0348 %0000110: 205b 6b48 7a5d 3a20 2020 2035 3030 300d   [kHz]:    5000.
0349 %0000120: 0a50 6572 696f 6473 3a20 2020 2020 2020  .Periods:
0350 %0000130: 2020 2020 2020 2020 2032 300d 0a46 7261           20..Fra
0351 %0000140: 6d65 733a 2020 2020 2020 2020 2020 2020  mes:
0352 %0000150: 2020 2020 2031 300d 0a0d 0a0d 0a
0353 %                                       8bc33e       10........>
0354 %0000160: 3f 0548633e bf20933d e192393d a568ea  ?.Hc>. .=..9=.h.
0355 %0000170: 3c 2530f63c 27e6043d 2043ad3c 25ce93  <%0.<'..= C.<%..
0356 %0000180: 3c aebcce3c 8714643d e3533d3e 65b6e1  <...<..d=.S=>e..
0357 %0000190: 3e 7a62103f 81c4143e 8c99813d 35921d  >zb.?...>...=5..
0358 %00001a0: 3d 8e6b0d3d 690cf93c 9910713c 3c9289  =.k.=i..<..q<<..
0359 %00001b0: 3c f6736f3c 2291453d ad1cab3d 386f15  <.so<".E=...=8o.
0360 %00001c0: 3e e82a143f 2a952e3f f568493e f08a8c  >.*.?*..?.hI>...
0361 %00001d0: 3d e43e0e3d 7040253d 19f4af3c 67fd93  =.>.=p@%=...<g..
0362 
0363 function vv = sheffield_readdata( fname );
0364    fid=fopen(fname,'rb','ieee-le');
0365    draw=fread(fid,[104,inf],'float32');
0366    fclose(fid);
0367    ldat = size(draw,2);
0368 
0369    [x,y]= meshgrid( 1:16, 1:16);
0370    idxm= y-x;
0371  % HW gain table
0372    gtbl = [0,849,213,87,45,28,21,19,21,28,45,87,213,849,0];
0373    idxm(idxm<=0)= 1;
0374    idxm= gtbl(idxm);
0375 
0376  % Multiply by gains
0377    draw = draw .* (idxm(idxm>0) * ones(1,ldat)); 
0378 
0379    vv= zeros(16*16, size(draw,2));
0380    vv(find(idxm),:) = draw;
0381 
0382  % Add reciprocal measurements
0383    vv= reshape(vv,[16,16,ldat]);
0384    vv= vv + permute( vv, [2,1,3]);
0385    vv= reshape(vv,[16*16,ldat]);
0386 
0387    
0388 
0389 function [vv,curr,volt,auxdata] = mceit_readdata( fname );
0390 
0391    fid= fopen(fname,'rb');
0392    d= fread(fid,inf,'float');
0393    fclose(fid);
0394 
0395    if rem( length(d), 256) ~=0
0396       eidors_msg('File length strange - cropping file',0);
0397       d=d(1:floor(length(d)/256)*256);
0398    end
0399 
0400    dd= reshape( d, 256, length(d)/256);
0401    if dd(39,1)==0 %104 measurements per frame
0402        dd=transfer104_208(dd);
0403    end
0404    vv= untwist(dd);
0405 
0406    curr=0.00512*dd(209:224,:);  % Amps
0407    volt=12*dd(225:240,:); %Vrms
0408    auxdata= dd(241:255,:);
0409    auxdata= auxdata(:);
0410    %input impedance=voltage./current-440;        Ohm
0411 
0412 function array208=transfer104_208(array104),
0413 % The order in the 208-vector is
0414 % Index   Inject    Measure Pair
0415 % 1         1-2        3-4        (alternate pair is 3-4, 1-2, at location 39)
0416 % 2         1-2        4-5
0417 % 3         1-2        5-6
0418 % ¡­
0419 % 13        1-2        15-16
0420 % 14        2-3        4-5
0421 % 15        2-3        5-6
0422 % ¡­
0423 % 26        2-3        16-1
0424 % 27        3-4        5-6
0425 % ¡­
0426 % 39        3-4        1-2
0427 % 40        4-5        6-7
0428 % ¡­
0429 
0430 ind=[39,51,63,75,87,99,111,123,135,147,159,171,183,52,64,76, ...
0431      88,100,112,124,136,148,160,172,184,196,65,77,89,101,113, ...
0432      125,137,149,161,173,185,197,78,90,102,114,126,138,150,162, ...
0433      174,186,198,91,103,115,127,139,151,163,175,187,199,104,116, ...
0434      128,140,152,164,176,188,200,117,129,141,153,165,177,189, ...
0435      201,130,142,154,166,178,190,202,143,155,167,179,191,203, ...
0436      156,168,180,192,204,169,181,193,205,182,194,206,195,207,208];
0437 ro=[1:13, 14:26, 27:38, 40:50, 53:62, 66:74, 79:86, 92:98, ...
0438     105:110,118:122,131:134,144:146,157:158,170:170];
0439 
0440 [x,y]=size(array104);
0441 if x~=256 && y~=256,
0442     eidors_msg(['eidors_readdata: expectingin an input array ', ...
0443                 'of size 208*n']);
0444     return;
0445 elseif y==256,
0446     array104=array104';
0447     y=x;
0448 end
0449 array208=array104;
0450 for i=1:y,
0451     array208(ind,i)=array104(ro,i);    
0452 end
0453    
0454 
0455 % measurements rotate with stimulation, we untwist them
0456 function vv= untwist(dd);
0457    elec=16;
0458    pos_i= [0,1];
0459    ELS= rem(rem(0:elec^2-1,elec) - ...
0460         floor((0:elec^2-1)/elec)+elec,elec)';
0461    ELS=~any(rem( elec+[-1 0 [-1 0]+pos_i*[-1;1] ] ,elec)' ...
0462         *ones(1,elec^2)==ones(4,1)*ELS')';
0463    twist= [               0+(1:13), ...
0464                          13+(1:13), ...
0465            39-(0:-1:0),  26+(1:12), ...
0466            52-(1:-1:0),  39+(1:11), ...
0467            65-(2:-1:0),  52+(1:10), ...
0468            78-(3:-1:0),  65+(1: 9), ...
0469            91-(4:-1:0),  78+(1: 8), ...
0470           104-(5:-1:0),  91+(1: 7), ...
0471           117-(6:-1:0), 104+(1: 6), ...
0472           130-(7:-1:0), 117+(1: 5), ...
0473           143-(8:-1:0), 130+(1: 4), ...
0474           156-(9:-1:0), 143+(1: 3), ...
0475           169-(10:-1:0),156+(1: 2), ...
0476           182-(11:-1:0),169+(1: 1), ...
0477           195-(12:-1:0), ...
0478           208-(12:-1:0) ];
0479     vv= zeros(256,size(dd,2));
0480     vv(ELS,:)= dd(twist,:);
0481    %vv= dd(1:208,:);
0482 
0483 % Read data from p2k files (I T S system)
0484 % FIXME: this code is very rough, it works for
0485 %   only eight ring data records
0486 function  vv = its_readdata( fname ) 
0487    fid= fopen( fname, 'rb', 'ieee-le');
0488    vv=[];
0489 
0490    % don't know how to interpret header
0491    header= fread(fid, 880, 'uchar');
0492    frameno= 0;
0493    rings= 8;
0494    while( ~feof(fid) )
0495        frameno= frameno+1;
0496        % don't know how to interpret frame header
0497        framehdr= fread(fid, 40);
0498        data= fread(fid, 104*rings, 'double');
0499        vv= [vv, data];
0500    end
0501 
0502    if 0 % convert a ring to 208
0503       ringno= 1;
0504       ld= size(vv,2);
0505       vx= [zeros(1,ld);vv( ringno*104 + (-103:0) ,: )];
0506       idx= ptr104_208;
0507       vv= vx(idx+1,:);
0508    end
0509 
0510 % pointer to convert from 104 to 208 meas patterns
0511 function idx= ptr104_208;
0512     idx= zeros(16);
0513     idx(1,:)= [0,0,1:13,0];
0514     ofs= 13;
0515     for i= 2:14
0516         mm= 15-i;
0517         idx(i,:) = [zeros(1,i+1), (1:mm)+ofs ];
0518         ofs= ofs+ mm;
0519     end
0520 
0521     idx= idx + idx';
0522     
0523 function vv = iirc_readdata( fname );
0524     fid= fopen( fname, 'r');
0525     while ~feof(fid)
0526        line = fgetl(fid);
0527        if isempty(line)
0528            continue;
0529        end
0530        
0531        num= regexp(line,'Channel : (\d+)');
0532        if ~isempty(num)
0533            channels= str2num( line(num(1):end ) );
0534            continue;
0535        end
0536        
0537        num= regexp(line,'Frequency : (\d+)kHz');
0538        if ~isempty(num)
0539            freqency= str2num( line(num(1):end ) );
0540            continue;
0541        end
0542 
0543        num= regexp(line,'Scan Method : (\w+)');
0544        if ~isempty(num)
0545            scan_method=  line(num(1):end );
0546            continue;
0547        end
0548 
0549        num= regexp(line,'Study : (\w+)');
0550        if ~isempty(num)
0551            study=  line(num(1):end);
0552            continue;
0553        end
0554            
0555        if strcmp(line,'Data');
0556            data= fscanf(fid,'%f',[4,inf])';
0557            continue;
0558        end
0559     end
0560     vv= data(:,1) + 1i*data(:,2);
0561     if length(vv) ~= channels^2
0562         error('eidors_readdata: data length wrong')
0563     end
0564 
0565 % stimulation is the fwd_model stimulation data structure
0566 % meas_select indicates if data is NOT measures on current electrode
0567 function [stimulations,meas_select] = UCT_sequence_file( fname );
0568    % (c) Tim Long
0569    % 21 January 2005
0570    % University of Cape Town
0571 
0572 
0573    % open the file
0574    fid = fopen(fname, 'rt');
0575 
0576    % check to see if file opened ok
0577    if fid == -1
0578          errordlg('File not found','ERROR')  
0579          return;
0580    end
0581 
0582 
0583    tline = fgetl(fid);             % get the spacer at top of text file
0584 
0585    % the measurement and injection pattern is stored as follows:
0586    % I1V1:  db #$00,#$0F,#$00,#$00,#$10,#$00,#$00,#$21,#$00 ...
0587    % I2V2:  db #$11,#$0F,#$11,#$11,#$10,#$11,#$11,#$21,#$11 ...
0588    % etc
0589    % need to put all the bytes in a vector
0590 
0591    % tokenlist will store the list of bytes as strings
0592    tokenlist = [];
0593 
0594 
0595    tline = fgetl(fid);             % get first line of data
0596 
0597    while length(tline) ~= 0
0598        
0599        % the first few characters in the line are junk
0600        rem = tline(11:end);            % extract only useful data
0601        
0602        % extract each byte
0603        while length(rem) ~=0
0604            [token, rem] = strtok(rem, ',');
0605            tokenlist = [tokenlist; token];
0606        end
0607        
0608        % get the next line in sequence file
0609        tline = fgetl(fid);
0610    end
0611 
0612    fclose(fid);
0613 
0614    % got everything in string form... need to covert to number format
0615 
0616    drive_lay = [];
0617    drive_elec = [];
0618    sense_lay = [];
0619 
0620    injection_no = 1;
0621    % for each injection
0622    for i=1:3:length(tokenlist)
0623        
0624        % get injection layer
0625        tsource_layer = tokenlist(i,3);
0626        tsink_layer = tokenlist(i,4);
0627        source_layer = sscanf(tsource_layer, '%x');
0628        sink_layer = sscanf(tsink_layer, '%x');
0629        
0630        drive_lay = [drive_lay; [source_layer sink_layer]];
0631          
0632        
0633        % get drive pair
0634        tsource_elec = tokenlist(i+1,3);
0635        tsink_elec = tokenlist(i+1,4);
0636        source_elec = sscanf(tsource_elec, '%x');
0637        sink_elec = sscanf(tsink_elec, '%x');
0638        
0639        drive_elec = [drive_elec; [source_elec sink_elec]];
0640        
0641        
0642        % get sense layer pair
0643        tpos_layer = tokenlist(i+2,3);
0644        tneg_layer = tokenlist(i+2,4);
0645        pos_sense_layer = sscanf(tpos_layer, '%x');
0646        neg_sense_layer = sscanf(tneg_layer, '%x');
0647        
0648        sense_lay = [sense_lay; [pos_sense_layer neg_sense_layer]];
0649    end
0650 
0651    n_elec = size(sense_lay,1);
0652    n_inj  = size(drive_lay,1);       % every injection
0653    elecs_per_plane = 16; % FIXED FOR THE UCT DEVICE
0654    raw_index = 0;
0655    meas_select = [];
0656 
0657    for i=1:n_inj      % for every injection
0658        stimulations(i).stimulation= 'mA';
0659        
0660        % find the injection electrodes
0661        e_inj_p = drive_lay(i, 1) * elecs_per_plane + drive_elec(i,1) + 1;
0662        e_inj_n = drive_lay(i, 2) * elecs_per_plane + drive_elec(i,2) + 1;
0663 
0664        % create the stimulation pattern for this injection
0665        inj = zeros(n_elec, 1);
0666        inj(e_inj_p) = 1;
0667        inj(e_inj_n) = -1;
0668        stimulations(i).stim_pattern = inj;
0669      
0670        % the UCT instrument always makes 16 measurements per injection.
0671        % the +ve and -ve electrodes are always adjacent, but might be on
0672        % different planes
0673        meas_pat = [];
0674        for e = 0:15
0675            raw_index = raw_index + 1;
0676            meas = zeros(1, n_elec);   % the measurement electrodes for this sample
0677 
0678            % find the measurement electrodes for this measurement (+ve elec is
0679            % next to -ve electrode)
0680            e_meas_p = sense_lay(i,1) * elecs_per_plane + mod(e+1,elecs_per_plane) + 1;
0681            e_meas_n = sense_lay(i,2) * elecs_per_plane + e + 1;
0682 
0683            % if either of the drive electrodes are equal to any of the sense
0684            % electrodes, we must not include this sample
0685            if any( e_meas_p == [e_inj_p, e_inj_n] ) | ...
0686               any( e_meas_n == [e_inj_p, e_inj_n] ) 
0687                continue;
0688            end
0689 
0690            meas(e_meas_p) = -1;
0691            meas(e_meas_n) = 1;
0692            meas_select = [meas_select;raw_index];
0693            
0694            % add this measurement to the measurement pattern
0695            meas_pat = [meas_pat; meas];
0696 
0697        end     % for each injection there are actually only 13-16 measurements
0698        stimulations(i).meas_pattern = sparse(meas_pat);
0699    end
0700 
0701 function [cur_data,no_cur_data] = UCT_calibration_file( fname );
0702    fid = fopen(fname, 'rb');
0703    mag_num = fread(fid, 1, 'int32');
0704    version = fread(fid, 1, 'int32');
0705 %  comments = fread(fid, 2048, 'char'); MATLAB CHANGED CHAR to 2 bytes
0706    comments = fread(fid, 2048, 'uint8');
0707    comments = setstr(comments');
0708    no_of_layers = fread(fid, 1, 'int32');
0709    uppa_dac = fread(fid, 1, 'float64');
0710    lowa_dac = fread(fid, 1, 'float64');
0711    raw_frame_size = fread(fid, 1, 'int32');
0712 
0713    no_cur_data = [];
0714    cur_data = [];
0715 
0716    for i=1:no_of_layers
0717        no_cur_data = [no_cur_data;fread(fid, raw_frame_size, 'float64')];
0718        cur_data = [cur_data;fread(fid, raw_frame_size, 'float64')];
0719    end
0720    fclose(fid);
0721 
0722 %  no_cur_data = UCT_ShuffleData( no_cur_data);
0723 %  cur_data =    UCT_ShuffleData( cur_data);
0724 
0725 function [v_raw] = UCT_LoadDataFrame(infilename)
0726 % [v_raw] = LoadDataFrame(infilename)
0727 %
0728 % Loads the data from the tomography file
0729 %
0730 % (c) Tim Long
0731 % 21 January 2005
0732 % University of Cape Town
0733 
0734 
0735 % Open the file
0736 fid = fopen(infilename, 'rb');
0737 
0738 if fid == -1
0739       errordlg('File not found','ERROR')  
0740       return;
0741 end
0742 
0743 %%% new file changes
0744 magic_number = fread(fid, 1, 'uint32');
0745 
0746 if magic_number ~= 2290649224
0747    disp('UCT File: wrong file type'); 
0748 end
0749 version = fread(fid, 1, 'uint32');
0750 foffset = fread(fid, 1, 'uint32');
0751 no_of_layers = fread(fid, 1, 'uint32');
0752 frame_size = fread(fid, 1, 'uint32');
0753 fseek(fid, foffset-8, 'cof');
0754 
0755 %%% end of new file changes
0756 %%% old file stuff
0757 % no_of_layers = fread(fid, 1, 'uint32');
0758 % frame_size = fread(fid, 1, 'uint32');
0759 
0760 frame_no = fread(fid, 1, 'uint32');
0761 
0762 v_raw = [];
0763 while feof(fid)==0
0764     v_raw = [v_raw, fread(fid, frame_size*no_of_layers, 'float64')];
0765     frame_no = fread(fid, 1, 'uint32');
0766 end
0767 
0768 fclose(fid);
0769 % UCT_ShuffleData??
0770 
0771 % Read data from the file format develped by Pascal
0772 % Gaggero at CSEM Landquart in Switzerland.
0773 function [vv] = landquart1_readdata( fname );
0774    FrameSize = 1024;
0775 
0776    enableCounter=1;
0777    counter=0;
0778 
0779    fid=fopen(fname);
0780 
0781    codec=fread(fid,1,'int'); %read codec
0782    if codec~=1; error('Codec unexpected value'); end
0783 
0784    nbFileInHeader=fread(fid,1,'int'); %read codec
0785    
0786    for i=1:nbFileInHeader
0787        lenghtFile = fread(fid,1,'int64'); 
0788        jnk= fread(fid,lenghtFile,'int8');% discard the file
0789    end
0790    
0791    
0792    vv= []; 
0793    while 1; 
0794        type=fread(fid,1,'int'); %check if not End of file
0795        if type == -1; break ; end
0796        
0797        nel= fread(fid,1,'int');
0798        sz= fread(fid,1,'int');
0799        iq=fread(fid,[2,sz/8],'int')';
0800 
0801        vv = [vv, iq*[1;1j]]; %I+ 1j*Q vectors
0802        
0803        for j=1:nel-1
0804            sz= fread(fid,1,'int');
0805            jnk = fread(fid,sz/4,'int');
0806        end
0807        
0808    end
0809    
0810 % Read data from the file format develped by Swisstom, Landquart, Switzerland.
0811 function [vv] = landquart2_readdata( fname )
0812    [fid msg]= fopen(fname,'r','ieee-be','UTF-8');
0813    try
0814       format_version = fread(fid,1,'int32','ieee-be');
0815       if format_version ~= 3
0816          error('unsupported file format version');
0817       else
0818          header_size = 2264; 
0819          fseek(fid,header_size + 8,'bof');
0820          %%% get frame size and length of payload at end of frame
0821          frame_length = fread(fid, 1, 'int32', 'ieee-be') + 12;
0822          %%% move back to start of 1. frame
0823          fseek(fid,header_size,'bof');
0824          
0825          %%% Read frames
0826          i = 1;
0827          while fseek(fid, 1,'cof') ~= -1
0828             fseek(fid, -1,'cof');
0829             % drop frame header and some payload:
0830             % (12 + 60) + 12 + 256 = 340
0831             fseek(fid,340, 'cof');
0832             iqPayload(:,i) = fread(fid,2048,'int32','ieee-le');
0833             fseek(fid,header_size + i*frame_length,'bof');
0834             i = i +1;
0835          end
0836 
0837       end
0838    catch err
0839       fclose(fid);
0840       rethrow(err);
0841    end
0842    fclose(fid);
0843 
0844    % this is just a simple guess
0845    amplitudeFactor = 2.048 / (2^20 * 360 * 1000);
0846    vv = amplitudeFactor * (iqPayload(1:2:end,:) + 1i*iqPayload(2:2:end,:));
0847 
0848 % Read data from the file format develped by Swisstom, Landquart, Switzerland.
0849 
0850 function [vv] = landquart4_readdata( fname )
0851    [fid msg]= fopen(fname,'r','ieee-le','UTF-8');
0852    try
0853       format_version = fread(fid,1,'int32','ieee-le');
0854       if format_version ~= 4
0855          error('unsupported file format version');
0856       else
0857          header_size = fread(fid,1,'int32', 'ieee-le'); 
0858          
0859          fseek(fid,header_size + 12,'bof');
0860          %%% get frame size and length of payload at end of frame
0861          frame_length = fread(fid, 1, 'int32', 'ieee-le') + 16;
0862          %%% move back to start of 1. frame
0863          fseek(fid,header_size,'bof');
0864          
0865          %%% Read frames
0866          i = 1;
0867          while fseek(fid, 1,'cof') ~= -1
0868             fseek(fid, -1,'cof');
0869             % drop frame header and some payload:
0870             % (16 + 60) + 12 + 256 = 344
0871             fseek(fid,344, 'cof');
0872             iqPayload(:,i) = fread(fid,2048,'int32','ieee-le');
0873             fseek(fid,header_size + i*frame_length,'bof');
0874             i = i +1;
0875          end
0876 
0877       end
0878    catch err
0879       fclose(fid);
0880       rethrow(err);
0881    end
0882    fclose(fid);
0883 
0884    % this is just a simple guess
0885    amplitudeFactor = 2.048 / (2^20 * 360 * 1000);
0886    vv = amplitudeFactor * (iqPayload(1:2:end,:) + 1i*iqPayload(2:2:end,:));
0887    
0888 %  Output: [encodepage] = eidors_readdata( path,'DIX_encode');
0889 %   where path= '/path/to/Criptografa_New.dll' (provided with system)
0890 function [vv] = dixtal_read_codepage( fname );
0891    %Fname must be the Criptografa_New dll
0892    fid = fopen(fname,'rb');
0893    b1234= fread(fid, [1,4], 'uint8');
0894    if b1234~= [77, 90, 144, 0];
0895       error('This does not appear to be the correct Dll');
0896    end
0897    fseek(fid, hex2dec('e00'), 'bof');
0898    encodepage1 = fread(fid, hex2dec('1000'), 'uint8');
0899    fclose(fid);
0900    encodepage1 = flipud(reshape(encodepage1,4,[]));
0901    vv = encodepage1(:);
0902 
0903 %        format = 'DIXTAL'
0904 %           - output: [vv] = eidors_readdata( fname, 'DIXTAL', encodepage );
0905 function [vv] = dixtal_read_data( file_name, frame_range, encodepage1 );
0906 
0907    % Get encodepage from the file
0908    fid=fopen(file_name,'rb');
0909    n1 = fgetl(fid);
0910    n2 = fgetl(fid);
0911    n3 = fgetl(fid); 
0912    nframes = str2num( n3 );
0913    
0914    fseek(fid, hex2dec('800'), 'bof');
0915    encodepage2 = fread(fid, hex2dec('1000'), 'uint8');
0916    fclose(fid);
0917 
0918    encodepage = bitxor(encodepage1, encodepage2);
0919 
0920    % Read the file
0921    dplen = hex2dec('1090');
0922    start = hex2dec('2800');
0923    fid=fopen(file_name,'rb');
0924    if isempty(frame_range)
0925       frame_range = 0:nframes;
0926    else
0927       frame_range = frame_range - 1;
0928       frame_range(frame_range>nframes) = [];
0929    end      
0930    vv= zeros(dplen/4, length(frame_range));
0931    k= 1;
0932    for i = frame_range
0933       status= fseek(fid, start + i*dplen, 'bof');
0934       if status~=0; break; end
0935       datapage = fread(fid, dplen, 'uint8');
0936       if length(datapage)== 0; 
0937          vv= vv(:,1:k-1); break; % frame number inside file is wrong
0938       end
0939       vv(:,k) = proc_dixtal_data( datapage, encodepage );
0940       k=k+1;
0941    end
0942    fclose(fid);
0943 
0944 
0945 function  data = proc_dixtal_data( b, encodepage );
0946 
0947    cryptolen = 1024*4;
0948    b(1:cryptolen) = bitxor(b(1:cryptolen), encodepage);
0949 
0950    b1 = b(1:4:end,:); %b1 = bitand(b1,127);
0951    b2 = b(2:4:end,:);
0952    b3 = b(3:4:end,:);
0953    b4 = b(4:4:end,:);
0954    bb  = [b1,b2,b3,b4]*(256.^[3;2;1;0]);
0955 
0956    sgnbit = bitand( bb,  2^31);
0957    sgnbit = +1*(sgnbit == 0) - 1*(sgnbit == 2^31);
0958 
0959    expbit = bitand( bb, 255*2^23 ) /2^23; 
0960    expbit = 2.^(expbit - 127);
0961 
0962    fracbt = bitand( bb, 2^23-1)/2^23 + 1;
0963    data = sgnbit .* expbit .* fracbt;
0964 
0965 function [fr] = read_draeger_header( filename );
0966 %READ_DRAEGER_HEADER   Opens and reads the header portion of the Draeger .eit
0967 %file. Current parameter returned is frame rate.
0968 %
0969 % function [fr,s] = ReadDraegerHeader(filename)
0970 % fr        double      scalar      frame rate in hertz
0971 
0972 % Determine file version
0973 K0 = 2048;   % bytes to read in char class
0974 K1='Framerate [Hz]:'; % Draeger frame rate line in header
0975 K2=15;                % Length of K1 string
0976 K3=13;                % Following F1 Field for delimiter (look for ^M)
0977 
0978 % Open file for reading in little-endian form
0979 fid = fopen(filename,'r','l');
0980 if fid == -1
0981     error('Error read_draeger_header: can''t open file');
0982 end
0983 header = fread(fid,K0,'*char')';
0984 index = strfind(header,K1);
0985 if ~isempty(index)
0986     [tok,rem]= strtok(header(index+K2:end),K3);
0987     fr = str2num(tok); 
0988 else
0989     error('Error read_draeger_header: frame rate unspecified');
0990 end
0991 
0992 if isempty(fr)
0993     error('Error read_draeger_header: Frame rate could not be read');
0994 end
0995 
0996 fclose(fid);
0997 
0998 
0999 function vd = read_draeger_file(filename)
1000 %READDRAEGERFILE   Opens and reads a Draeger .eit file. Returns an array
1001 %containing the voltage data arranged per frame.
1002 %
1003 % Input:
1004 % filename  char        1xN
1005 %
1006 % Output:
1007 % vd        double      MxN         M = volt measurements per frame (mV)
1008 %                                   N = number of frames
1009 
1010 
1011    ff = dir(filename);
1012    filelen = ff.bytes;
1013 
1014    % Open file for reading in little-endian form
1015    fid = fopen(filename,'r','l');
1016    if fid == -1
1017        error('Error read_draeger_file: file could not be opened');
1018    end
1019 
1020    % Determine file version
1021    K0 = 128;   % bytes to read in char class
1022 
1023    header = fread(fid,K0,'*char')';
1024 if 0 % THis analysis doesn't seem useful
1025    if ~isempty(strfind(header,'V3.2'))
1026        version = '3.2';
1027    elseif ~isempty(strfind(header,'V4.01'))
1028        version = '4.01'; 
1029    elseif ~isempty(strfind(header,'V1.10')) % 2014!!
1030        version = '1.10'; 
1031    else
1032        error('Error read_draeger_file: unknown file version');
1033    end
1034 end
1035 
1036    % Read Format version
1037    fseek(fid,0,-1); % beginning of file
1038    hdr = fread(fid, 8, 'uint8');
1039    offset1 =  (256.^(0:3))*hdr(5:8) + 16; % "good" data starts at offset #2
1040    switch sprintf('%02X-',header(1:4));
1041       case '20-00-00-00-';
1042          type='Draeger v32';  SPC = 5200;
1043       case '33-00-00-00-';
1044          type='Draeger v51';  SPC = 5495;
1045       otherwise;
1046          error('File "%s" is format version %d, which we don''t know how to read', ...
1047                hdr(1))
1048    end
1049 
1050    len = floor( (filelen-offset1)/SPC ) - 1;
1051    vd= zeros(512,len);
1052    ss= offset1 + (0:len)*SPC;
1053 
1054    for k= 1:length(ss);
1055       if fseek(fid,ss(k),-1)<0; break; end
1056       vd(:,k)=fread(fid,512,'double', 0, 'ieee-le');
1057    end
1058 
1059    % End of function
1060    fclose(fid);
1061

Generated on Tue 12-May-2015 16:00:42 by m2html © 2005