find_electrode_bdy

 FIND_ELECTRODE_BDY: find the boundary index area for electrode
 Input:
   bdy => boundary (from fwd_model.boundary) bdy simplices x nodes index 
   vtx => node pts (from fwd_model.nodes)
   elec_nodes => index of nodes in the electrode
 Output:
   bdy_idx  => vector of boundary simplices in this electrode
   bdy_area => boundary area of each simplex in bdy_idx
  if the nodes in the electrode are points, then
   bdy_area is the area corresponding to these point electrodes

 if [bdy_idx, bdy_area] = find_electrode_bdy( bdy, vtx, [])
   calculate bdy_idx=[] and bdy_area for all on bdy

0001 function [bdy_idx, bdy_area] = find_electrode_bdy( bdy, vtx, elec_nodes)
0002 % FIND_ELECTRODE_BDY: find the boundary index area for electrode
0003 % Input:
0004 %   bdy => boundary (from fwd_model.boundary) bdy simplices x nodes index
0005 %   vtx => node pts (from fwd_model.nodes)
0006 %   elec_nodes => index of nodes in the electrode
0007 % Output:
0008 %   bdy_idx  => vector of boundary simplices in this electrode
0009 %   bdy_area => boundary area of each simplex in bdy_idx
0010 %  if the nodes in the electrode are points, then
0011 %   bdy_area is the area corresponding to these point electrodes
0012 %
0013 % if [bdy_idx, bdy_area] = find_electrode_bdy( bdy, vtx, [])
0014 %   calculate bdy_idx=[] and bdy_area for all on bdy
0015 
0016 % (C) 2008 Andy Adler. License: GPL version 2 or version 3
0017 % $Id: find_electrode_bdy.m 5960 2019-06-05 02:24:24Z aadler $
0018 
0019 % special processing if elec_nodes==[]
0020 if isempty(elec_nodes);
0021    bdy_idx = [];
0022    l_bdy= size(bdy,1);
0023    bdy_area = zeros(l_bdy,1);
0024    for i=1:l_bdy
0025       bdy_area(i) = tria_area( vtx(bdy(i,:),:) );
0026    end
0027    return;  
0028 end
0029 
0030 [bdy_idx, point] = find_bdy_idx( bdy, elec_nodes);
0031 if nargout==1; return;end
0032 
0033 l_bdy_idx= length(bdy_idx);
0034 l_point= length(point);
0035 
0036 if l_bdy_idx > 0 && l_point ==0
0037    bdy_area= zeros(1, l_bdy_idx);
0038 
0039    for i=1:l_bdy_idx
0040       bdy_nodes   = bdy(bdy_idx(i),:);
0041       bdy_area(i) = tria_area( vtx(bdy_nodes,:) );
0042    end
0043 elseif l_bdy_idx == 0 && l_point >0
0044    dims = size(bdy,2); 
0045    l_point= length(point);
0046    bdy_area = zeros(l_point,1);
0047    for i=1:l_point
0048      ff= find( any(bdy== point(i),2) );
0049      this_area= 0;
0050      for ffp=ff(:)'
0051         xyz= vtx( bdy(ffp,:),:);
0052         this_area= this_area + tria_area( xyz );
0053      end
0054      bdy_area(i)= bdy_area(i) + this_area/dims;
0055    end
0056 else
0057 keyboard
0058    error(['Can''t model this electrode. ' ...
0059           'It has %d CEM and %d point nodes (e.g. nodes %d, %d)'], ...
0060             l_bdy_idx, l_point, bdy_idx(1), point(1) )
0061 end
0062    
0063 
0064 
0065 % Find which boundary nodes are part of a
0066 % complete simplex. These nodes can be made
0067 % into complete electrode model elements
0068 %
0069 % Any nodes left will be treated as point electrodes
0070 function [ffb,unused] = find_bdy_idx( bdy, elec_nodes);
0071    bdy_els = zeros(size(bdy,1),1);
0072    elec_nodes = unique(elec_nodes);
0073    for nd= elec_nodes(:)'
0074       bdy_els = bdy_els + any(bdy==nd,2);
0075    end
0076    ffb = find(bdy_els == size(bdy,2));
0077 
0078 %  find if all nodes are used
0079    used_nodes= reshape( bdy(ffb,:), 1,[]);
0080    unused=     setdiff( elec_nodes, used_nodes);
0081    
0082 
0083 % bdy points is [x1,y1,z1;x2,y2,z2; etc]
0084 function area= tria_area( bdy_pts ); 
0085    vectors= diff(bdy_pts); 
0086    if size(vectors,1)==2
0087       vectors= cross(vectors(1,:),vectors(2,:));
0088    end
0089    d= size(bdy_pts,1);
0090    area= sqrt( sum(vectors.^2) )/( d-1 );
0091