mk_fmdl_from_nodes

 MK_FMDL_FROM_NODES: create fmdl from nodes
 fwd_mdl= mk_fmdl_from_nodes( vtx, elec_nodes, z_contact, name)

 Create a fwd_model by delaynay triangularization of vtx nodes
  vtx:         matrix of nodes in model (can (or not) include elec_nodes)
  elec_nodes:  cell of matrix N x [x,y,{z}] for each electrode
  z_contact:   vector or scalar electrode contact impedance
  name:        name for eidors object

  fwd_mdl:     eidors format fwd_model

 Limitations: assumes a convex model

0001 function mdl= mk_fmdl_from_nodes( vtx, elec_nodes, z_contact, name)
0002 % MK_FMDL_FROM_NODES: create fmdl from nodes
0003 % fwd_mdl= mk_fmdl_from_nodes( vtx, elec_nodes, z_contact, name)
0004 %
0005 % Create a fwd_model by delaynay triangularization of vtx nodes
0006 %  vtx:         matrix of nodes in model (can (or not) include elec_nodes)
0007 %  elec_nodes:  cell of matrix N x [x,y,{z}] for each electrode
0008 %  z_contact:   vector or scalar electrode contact impedance
0009 %  name:        name for eidors object
0010 %
0011 %  fwd_mdl:     eidors format fwd_model
0012 %
0013 % Limitations: assumes a convex model
0014 
0015 % (C) 2008 Andy Adler. License: GPL version 2 or version 3
0016 % $Id: mk_fmdl_from_nodes.m 3600 2012-10-07 00:44:58Z aadler $
0017 
0018 mdl= eidors_obj('fwd_model', name);
0019 
0020 for i= 1:prod(size(elec_nodes)) 
0021    vtx= [vtx; elec_nodes{i}];
0022 end
0023 
0024 % Sort and remove duplicates
0025 vtx = unique(vtx, 'rows');
0026 
0027 %vtx=  perturb_vtx(vtx) ;
0028 %simp= delaunayn( vtx  );
0029 simp= delaunayn( perturb_vtx(vtx)  );
0030 % remove zero area elements
0031 keep= ones(size(simp,1),1);
0032 oo= ones(size(simp,2),1);
0033 for i=1:size(simp,1);
0034    area= abs(det([oo,vtx(simp(i,:),:)]));
0035    if area<1e-6; keep(i)=0;end
0036 end
0037 simp= simp(find(keep),:);
0038 
0039 % reorder elements so they're sorted by lowest node number
0040 [jnk,idx] = sort(mean(simp,2));
0041 simp= simp(idx,:);
0042 
0043 
0044 mdl.nodes    = vtx;
0045 mdl.elems    = simp;
0046 
0047 
0048 mdl.boundary = find_boundary( simp );
0049 mdl.gnd_node = 1;
0050 
0051 % Electrodes and z_contact
0052 
0053 n_elec= prod(size(elec_nodes));
0054 z_contact= z_contact.*ones(n_elec,1);
0055 for i= 1:n_elec
0056    [jnk, idxa]            = intersect( vtx, elec_nodes{i},'rows');
0057    electrodes(i).nodes    = idxa(:)';
0058    electrodes(i).z_contact= z_contact(i);
0059 end
0060 
0061 
0062 mdl.electrode =     electrodes;
0063 mdl.solve=          'eidors_default';
0064 mdl.jacobian=       'eidors_default';
0065 mdl.system_mat=     'eidors_default';
0066 mdl.normalize_measurements = 0;
0067 
0068 
0069 function vtx_perturb= perturb_vtx( vtx );
0070 if 0
0071    ctr= mean(vtx,1);
0072    v_ctr= vtx - ones(size(vtx,1),1) * ctr;
0073    r2 = (sum(v_ctr.^2,2)) .^ (1/4);
0074    vtx_perturb = vtx + .1* v_ctr ./ (r2*ones(1,size(vtx,2)));
0075 else
0076    scl= std(vtx,1);
0077    vtx_perturb = vtx + 1e-5*randn(size(vtx)).*(ones(size(vtx,1),1)*scl);
0078 end
0079 
0080    
0081 function vtx_perturb= perturb_vtx_try1( vtx );
0082    max_vtx = max(vtx);
0083    min_vtx = min(vtx);
0084    vtx_perturb = vtx + 0.05*randn(size(vtx));
0085    for d= 1:size(vtx,2)
0086       idx= (vtx(:,d) == min_vtx(d)) | ...
0087            (vtx(:,d) == max_vtx(d));
0088       vtx_perturb(idx,d) = vtx(idx,d); 
0089    end