meas_icov_rm_elecs

 MEAS_ICOV_RM_ELECS: remove electrodes from consideration
 meas_icov = meas_icov_rm_elecs( inv_model, elec_list )

 PARAMETERS:
 - elec_list = numbers of elecs to remove
     or as imdl.meas_icov_rm_elecs.elec_list

 - imdl.meas_icov_rm_elecs.exponent - exponent
 - imdl.meas_icov_rm_elecs.SNR      - SNR to add (default inf)
 - imdl.meas_icov_rm_elecs.replace_value (default 0)
      Default is to modify the current meas_icov value, if replace_value==1,
      then a new value is calculated without reference to the current

 meas_icov_rm_elecs can also accept a fwd_model parameter
 
 Reference Accounting for erroneous electrode data in EIT
 A. Adler Physiological Measurement, 25(1):227-238, 2004.

0001 function meas_icov = meas_icov_rm_elecs( imdl, elec_list)
0002 % MEAS_ICOV_RM_ELECS: remove electrodes from consideration
0003 % meas_icov = meas_icov_rm_elecs( inv_model, elec_list )
0004 %
0005 % PARAMETERS:
0006 % - elec_list = numbers of elecs to remove
0007 %     or as imdl.meas_icov_rm_elecs.elec_list
0008 %
0009 % - imdl.meas_icov_rm_elecs.exponent - exponent
0010 % - imdl.meas_icov_rm_elecs.SNR      - SNR to add (default inf)
0011 % - imdl.meas_icov_rm_elecs.replace_value (default 0)
0012 %      Default is to modify the current meas_icov value, if replace_value==1,
0013 %      then a new value is calculated without reference to the current
0014 %
0015 % meas_icov_rm_elecs can also accept a fwd_model parameter
0016 %
0017 % Reference Accounting for erroneous electrode data in EIT
0018 % A. Adler Physiological Measurement, 25(1):227-238, 2004.
0019 
0020 % (C) 2008 Andy Adler. License: GPL version 2 or version 3
0021 % $Id: meas_icov_rm_elecs.m 5112 2015-06-14 13:00:41Z aadler $
0022 
0023 if ischar(imdl) && strcmp(imdl,'UNIT_TEST'); do_unit_test; return; end
0024 
0025 switch imdl.type,
0026   case 'inv_model'; fmdl = imdl.fwd_model;
0027   case 'fwd_model'; fmdl = imdl;
0028                     imdl.meas_icov_rm_elecs.replace_value = 1;
0029   otherwise;        error('meas_icov_rm_elecs: require inv- or fwd-model');
0030 end
0031 
0032 if nargin< 2
0033    elec_list = imdl.meas_icov_rm_elecs.elec_list;
0034 end
0035 
0036      NSR = 0;
0037 try; NSR = (imdl.meas_icov_rm_elecs.SNR)^(-1);
0038 end
0039 
0040     exponent = 1;
0041 try;exponent = imdl.meas_icov_rm_elecs.exponent;
0042 end
0043 
0044     replace_value = 0;
0045 try;replace_value = imdl.meas_icov_rm_elecs.replace_value;
0046 end
0047  
0048 
0049 meas_icov = [];
0050 for stim = fmdl.stimulation(:)'
0051    mp = stim.meas_pattern;
0052    sp = stim.stim_pattern;
0053    icovi = ones(size(mp,1),1);
0054    if any(sp(elec_list) ~= 0)
0055       icovi = 0*icovi;
0056    else
0057       icovi = ~any( mp(:,elec_list) ~= 0, 2);
0058    end
0059 
0060    meas_icov = [meas_icov; icovi];
0061 end
0062 
0063 n = length(meas_icov);
0064 
0065 meas_icov(meas_icov == 0) = NSR;
0066 meas_icov = spdiags( meas_icov.^exponent, 0, n,n );
0067 
0068 if replace_value == 0
0069    meas_icov = calc_meas_icov(imdl)*meas_icov;
0070 end
0071 
0072 function do_unit_test
0073    imdl = mk_common_model('a2c0',4);
0074    covar = meas_icov_rm_elecs(imdl,1);
0075    unit_test_cmp('ire # 1', covar, zeros(size(covar))); 
0076 
0077    imdl = mk_common_model('a2c0',8);
0078    covar = meas_icov_rm_elecs(imdl,1);
0079    ff =    find( diag(covar)~= 1);
0080    ffcmp = [1;2;3;4;5;10;11;15;16;20;21;25;26;30;31;36;37;38;39;40];
0081    unit_test_cmp('ire # 2', ff,ffcmp);
0082 
0083    covar = meas_icov_rm_elecs(imdl.fwd_model,1);
0084    unit_test_cmp('ire # 3', ff,ffcmp);
0085 
0086    imdl.meas_icov_rm_elecs.elec_list = 1;
0087    covar = meas_icov_rm_elecs(imdl);
0088    unit_test_cmp('ire # 4', ff,ffcmp);
0089    ff =    find( diag(covar)~= 1);
0090 
0091    
0092    covar = meas_icov_rm_elecs(imdl,[1,2]);
0093    ff =    find( diag(covar)==1);
0094    ffcmp = [12;13;14;18;19;23;24;28;29;33;34;35];
0095    unit_test_cmp('ire # 5', ff,ffcmp);
0096 
0097    covar = meas_icov_rm_elecs(imdl,[1,2,3]);
0098    ff =    find( diag(covar)==1);
0099    ffcmp = [18;19;24;29;34;35];
0100    unit_test_cmp('ire # 6', ff,ffcmp);
0101 
0102    covar = meas_icov_rm_elecs(imdl,[1,2,3,4]);
0103    ff =    find( diag(covar)==1);
0104    ffcmp = [24;35];
0105    unit_test_cmp('ire # 7', ff,ffcmp);
0106 
0107 % NOW CHECK OTHER MODES
0108    imdl.meas_icov_rm_elecs.elec_list = 1:4;
0109    imdl.meas_icov_rm_elecs.exponent = 1;
0110    imdl.meas_icov_rm_elecs.SNR      = 0;
0111    covar=  diag( meas_icov_rm_elecs(imdl) );
0112    ffcmp = [24;35];
0113    unit_test_cmp('ire # 8', find(covar==1),ffcmp);
0114 
0115    imdl.meas_icov_rm_elecs.SNR      = 100;
0116    covar=  diag( meas_icov_rm_elecs(imdl) );
0117    unit_test_cmp('ire # 9', find(covar==1),ffcmp);
0118    ff =    find( covar~=1);
0119    unit_test_cmp('ire #10', covar(ff),1/100);
0120 
0121    imdl.meas_icov_rm_elecs.exponent = -1;
0122    covar=  diag( meas_icov_rm_elecs(imdl) );
0123    unit_test_cmp('ire #11', find(covar==1),ffcmp);
0124    ff =    find( covar~=1);
0125    unit_test_cmp('ire #12', covar(ff),100);
0126 
0127    imdl = mk_common_model('a2c0',8);
0128    imdl.meas_icov = spdiag((1:40)');
0129    covar=  diag( meas_icov_rm_elecs(imdl,[]) );
0130    unit_test_cmp('ire #13', covar,(1:40)');
0131 
0132    imdl.meas_icov_rm_elecs.replace_value = 0;
0133    covar=  diag( meas_icov_rm_elecs(imdl,[]) );
0134    unit_test_cmp('ire #14', covar,(1:40)');
0135 
0136    imdl.meas_icov_rm_elecs.replace_value = 1;
0137    covar=  diag( meas_icov_rm_elecs(imdl,[]) );
0138    unit_test_cmp('ire #15', covar,1);
0139