np_fwd_parameters

 NP_FWD_PARAMETERS: data= np_fwd_solve( fwd_model )
 Extract parameters from a 'fwd_model' struct which are 
 appropriate for Nick Polydorides EIDORS3D code
   param.n_elem   => number of elements
   param.n_elec   => number of electrodes
   param.n_node   => number of nodes (vertices)
   param.n_stim   => number of current stimulation patterns
   param.n_meas   => number of measurements (total)
   param.vtx      => vertex matrix
   param.simp     => connection matrix
   param.srf      => boundary triangles
   param.df       => vector of measurements for each current pattern
   param.elec     => nodes attached to each electrode
   param.zc       => vector of contact impedances
   param.indH     => electrodes used for each measurement
   param.I        => RHS (current term) for FEM solution
   param.Ib       => Current for electrodes
   param.perm_sym => 'sym' parameter
   param.gnd_ind  => node attached to ground
   param.normalize  => difference measurements normalized?

0001 function param = np_fwd_parameters( fwd_model )
0002 % NP_FWD_PARAMETERS: data= np_fwd_solve( fwd_model )
0003 % Extract parameters from a 'fwd_model' struct which are
0004 % appropriate for Nick Polydorides EIDORS3D code
0005 %   param.n_elem   => number of elements
0006 %   param.n_elec   => number of electrodes
0007 %   param.n_node   => number of nodes (vertices)
0008 %   param.n_stim   => number of current stimulation patterns
0009 %   param.n_meas   => number of measurements (total)
0010 %   param.vtx      => vertex matrix
0011 %   param.simp     => connection matrix
0012 %   param.srf      => boundary triangles
0013 %   param.df       => vector of measurements for each current pattern
0014 %   param.elec     => nodes attached to each electrode
0015 %   param.zc       => vector of contact impedances
0016 %   param.indH     => electrodes used for each measurement
0017 %   param.I        => RHS (current term) for FEM solution
0018 %   param.Ib       => Current for electrodes
0019 %   param.perm_sym => 'sym' parameter
0020 %   param.gnd_ind  => node attached to ground
0021 %   param.normalize  => difference measurements normalized?
0022 
0023 
0024 % (C) 2005 Andy Adler. License: GPL version 2 or version 3
0025 % $Id: np_fwd_parameters.m 5394 2017-04-12 15:10:30Z aadler $
0026 
0027 warning('EIDORS:deprecated','NP_FWD_PARAMETERS is deprecated as of 07-Jun-2012. ');
0028 
0029 try; switch fwd_model.type
0030         case 'fwd_model'; % do nothing
0031         case 'inv_model'; fwd_model = fwd_model.fwd_model;
0032         case 'image';     fwd_model = fwd_model.fwd_model;
0033         otherwise; error('np_fwd_parameters: requires fwd_model, inv_model or image object');
0034         end
0035 catch
0036      error('np_fwd_parameters: requires EIDORS object');
0037 end
0038 
0039 
0040 param = eidors_obj('get-cache', fwd_model, 'np_fwd_parameters');
0041 
0042 if ~isempty(param)
0043    eidors_msg('np_fwd_parameters: using cached value', 3);
0044    return
0045 end
0046 
0047 param = calc_param( fwd_model );
0048 
0049 eidors_obj('set-cache', fwd_model, 'np_fwd_parameters', param);
0050 eidors_msg('np_fwd_parameters: setting cached value', 3);
0051 
0052 % perform actual parameter calculation
0053 function param= calc_param( fwd_model );
0054 
0055 vtx= fwd_model.nodes;
0056 simp= fwd_model.elems;
0057 % calc num electrodes, nodes, stim_patterns
0058 n_elem= size(simp,1);
0059 n_elec=  length(fwd_model.electrode );
0060 n_node = size(fwd_model.nodes,1);
0061 n_stim = length(fwd_model.stimulation );
0062 n_meas = 0;
0063 
0064 % Recreate 'df' from fwd_model.stimulation
0065 df= zeros(n_stim,1);
0066 for i=1:n_stim;
0067     % Fixed by Steffen Kaufmann 06.03.2014 - should now work with non full
0068     % rank pattern
0069     % df(i) = size(fwd_model.stimulation(i).meas_pattern ,1);
0070     df(i) = sum(sum(abs(full(fwd_model.stimulation(i).meas_pattern)))) / 2;
0071     n_meas = n_meas + df(i);
0072 end
0073 
0074 elec= [];
0075 zc  = zeros(n_elec, 1);
0076 
0077 if isfield(fwd_model,'boundary')
0078     srf = fwd_model.boundary;
0079 else
0080     srf= find_boundary(simp);
0081 end
0082 
0083 max_elec_nodes=0;
0084 % get electrode parameters
0085 for i=1:n_elec
0086     elec_nodes= fwd_model.electrode(i).nodes;
0087     if length(elec_nodes)>1
0088        e_bdy  = bdy_with_nodes(srf,  elec_nodes );
0089        n_bdy  = srf(e_bdy,:)';
0090     else
0091        n_bdy= elec_nodes;
0092     end
0093 % elec is a series of nodes matching bdy faces
0094     en_list{i}= n_bdy(:)';
0095     if length(n_bdy) > max_elec_nodes
0096         max_elec_nodes = length(n_bdy);
0097     end
0098 
0099 % contact impedance
0100     zc(i)    = fwd_model.electrode(i).z_contact;
0101 end
0102 
0103 elec= zeros(n_elec, max_elec_nodes);
0104 for i=1:n_elec
0105     en= en_list{i};
0106     elec(i,1:length(en)) = en;
0107 end
0108 
0109 % Recreate 'indH' from fwd_model.stimulation
0110 indH= zeros(n_stim, 2);
0111 idx=0;
0112 for i=1:n_stim
0113    meas_pat= fwd_model.stimulation(i).meas_pattern';
0114 
0115    sourcepos= find(meas_pat(:)== 1);
0116    sourcepos= rem( sourcepos-1 , n_elec) + 1;
0117 
0118    sinkpos  = find(meas_pat(:)==-1);
0119    sinkpos  = rem( sinkpos  -1 , n_elec) + 1;
0120 
0121    indH( idx+(1:df(i)) , : ) = [sourcepos, sinkpos];
0122    idx= idx+ df(i);
0123 end
0124 
0125 % calculate FEM RHS matrix, i.e., the current patterns padded with zeroes
0126 I = zeros( n_elec + n_node, n_stim );
0127 idx=0;
0128 for i=1:n_stim
0129    I( n_node + (1:n_elec), i ) = ...
0130          fwd_model.stimulation(i).stim_pattern;
0131 end
0132 I(fwd_model.gnd_node,:) = 0;
0133 Ib= I( n_node + (1:n_elec), : );
0134 
0135 % pack into a parameter return list
0136 param.n_elem   = n_elem;
0137 param.n_elec   = n_elec;
0138 param.n_node   = n_node;
0139 param.n_stim   = n_stim;
0140 param.n_meas   = n_meas;
0141 param.vtx      = vtx;
0142 param.simp     = simp;
0143 param.srf      = srf;
0144 param.df       = df;
0145 param.elec     = elec;
0146 param.zc       = zc;
0147 param.indH     = indH;
0148 param.I        = I;
0149 param.Ib       = Ib;
0150 try
0151    param.perm_sym = fwd_model.np_fwd_solve.perm_sym;
0152 catch
0153    param.perm_sym = '{n}';
0154 end
0155 param.gnd_ind  = fwd_model.gnd_node;
0156 
0157 if isfield(fwd_model,'normalize_measurements')
0158    param.normalize = fwd_model.normalize_measurements;
0159 else
0160    param.normalize = 0;
0161 end
0162 
0163 
0164 
0165 % get boundary faces which match nodes
0166 function e_bdy  = bdy_with_nodes(bdy,  elec_nodes );
0167    mbdy= zeros(size(bdy));
0168    for n= elec_nodes(:)'
0169       mbdy= mbdy + (bdy == n); 
0170    end 
0171    e_bdy = find( all(mbdy') );
0172 
0173 % get boundary faces which match any node
0174 % Use this for point electrodes where there are no bdy faces
0175 % This is sort of an abuse of the model, but at least it can
0176 % produce a reasonable result for pt electrode mdls.
0177 
0178    if isempty(e_bdy)
0179       e_bdy = find( sum(mbdy')>=2 );
0180    end
0181    if isempty(e_bdy)
0182       e_bdy = find( any(mbdy') );
0183    end