join_models

 JOIN_MODELS: Join two fmdl structures to create one

 [fmdlu]= join_models(fmdl1, fmdl2, tol)
 fmdlu is the union of fmdl1 and fmdl2 with common nodes removed.
 tol is the minimum distance between nodes to consider them unique
 tol (default 1e-6) if not provided. 

 join_models only works on the elems, nodes and electrode fields;
 any c2f, mat_idx, stimulation structures must be redefined

 See also: MERGE_MESHES

0001 function [fmdlo]= join_models(fmdl1, fmdl2, tol)
0002 % JOIN_MODELS: Join two fmdl structures to create one
0003 %
0004 % [fmdlu]= join_models(fmdl1, fmdl2, tol)
0005 % fmdlu is the union of fmdl1 and fmdl2 with common nodes removed.
0006 % tol is the minimum distance between nodes to consider them unique
0007 % tol (default 1e-6) if not provided.
0008 %
0009 % join_models only works on the elems, nodes and electrode fields;
0010 % any c2f, mat_idx, stimulation structures must be redefined
0011 %
0012 % See also: MERGE_MESHES
0013 %
0014 
0015 % (C) 2018 Andy Adler. License: GPL version 2 or version 3
0016 % $Id: join_models.m 5848 2018-11-10 18:22:43Z aadler $
0017 
0018 if ischar(fmdl1) && strcmp(fmdl1,'UNIT_TEST'); do_unit_test; return; end
0019 
0020 if nargin==2; tol = 1e-5; end
0021 
0022 fmdlo = do_join_models(fmdl1,fmdl2,tol);
0023 
0024 function fmdlo = do_join_models(fmdl1,fmdl2,thresh);
0025    n1 = fmdl1.nodes;
0026    n2 = fmdl2.nodes;
0027    Dnodes = 0;
0028    nD = size(n1,2); %dimension
0029    for i=1:nD
0030       Dnodes = Dnodes + abs(bsxfun(@minus, n1(:,i), n2(:,i)'));
0031    end
0032    % idx of nodes in fmdl2 equal to nodes in fmdl1
0033    idx = find(any( Dnodes<thresh,1 ));
0034 
0035    nn1 = num_nodes(fmdl1);
0036    fmdlo = eidors_obj('fwd_model','joined model');
0037    fmdlo.nodes = [fmdl1.nodes;fmdl2.nodes];
0038    fmdlo.elems = [fmdl1.elems;fmdl2.elems+nn1];
0039    fmdlo.gnd_node = fmdl1.gnd_node;
0040    oidx = 1:num_nodes(fmdlo);
0041    nno = num_nodes(fmdlo);
0042    oidx(idx + nn1)= [];
0043    nidx = zeros(num_nodes(fmdlo),1);
0044    nidx(oidx) = 1:length(oidx);
0045    for i= idx(:)'
0046       ff = find(Dnodes(:,i)<thresh);
0047       if length(ff)~=1; error('Degenerate models with equal nodes'); end
0048       nidx(i+nn1) = ff;
0049    end
0050    fmdlo.nodes = fmdlo.nodes(oidx,:);
0051    fmdlo.elems = reshape(nidx(fmdlo.elems(:)),[],nD+1);
0052    if isfield(fmdl1,'electrode')
0053       fmdlo.electrode = fmdl1.electrode;
0054    else
0055       % SUPID MATLAB FORCES THIS SILLY INITALIZATION
0056       fmdlo.electrode = struct('nodes',{},'z_contact',{});
0057    end
0058    if isfield(fmdl2,'electrode')
0059       for i=1:length(fmdl2.electrode);
0060          eli = fmdl2.electrode(i);
0061          eli.nodes = nidx(eli.nodes + nn1);
0062          fmdlo.electrode(end+1) = eli;
0063       end
0064    end
0065    fmdlo.boundary = find_boundary(fmdlo.elems);
0066 
0067    for i=1:num_elecs(fmdlo)
0068       fmdlo.electrode(i).nodes = fmdlo.electrode(i).nodes(:)'; % force row
0069    end
0070 
0071 
0072 function do_unit_test
0073 
0074    subplot(221); [fmdl1,fmdl2]=do_unit_test_2D_mdls;
0075 
0076    fmdlo= join_models(fmdl1, fmdl2);
0077    do_testing('join_models-2Da',fmdlo,fmdl1,fmdl2)
0078    subplot(222); show_fem(fmdlo,[0,1,1]); axis off
0079 
0080    fmdl1 = rmfield(fmdl1,'electrode');
0081    fmdlo= join_models(fmdl1, fmdl2);
0082    do_testing('join_models-2Db',fmdlo,fmdl1,fmdl2)
0083 
0084    fmdl2 = rmfield(fmdl2,'electrode');
0085    fmdlo= join_models(fmdl1, fmdl2);
0086    do_testing('join_models-2Dc',fmdlo,fmdl1,fmdl2)
0087 
0088    subplot(223); [fmdl1,fmdl2]=do_unit_test_3D_mdls; view(0,63);
0089    fmdlo= join_models(fmdl1, fmdl2);
0090    do_testing('join_models-3Da',fmdlo,fmdl1,fmdl2)
0091    subplot(224); show_fem(fmdlo,[0,1,0]); axis off; view(0,63);
0092 
0093 
0094 function do_testing(txt,fmdlo,fmdl1,fmdl2)
0095    unit_test_cmp([txt,'-01'], ...
0096        size(fmdl1.elems,1) + size(fmdl2.elems,1), ...
0097        size(fmdlo.elems,1) );
0098    unit_test_cmp([txt,'-02'], ...
0099        size(fmdl1.nodes,1) + size(fmdl2.nodes,1) >= ...
0100        size(fmdlo.nodes,1), true );
0101 
0102    if ~isfield(fmdl1,'electrode'); fmdl1.electrode = struct([]); end
0103    if ~isfield(fmdl2,'electrode'); fmdl2.electrode = struct([]); end
0104    unit_test_cmp([txt,'-02'], ...
0105        length(fmdl1.electrode) + length(fmdl2.electrode), ...
0106        length(fmdlo.electrode));
0107 
0108 function [fmdl1,fmdl2]=do_unit_test_2D_mdls;
0109    imdl = mk_common_model('a2c0',8); fmdl= imdl.fwd_model;
0110    fmdl = crop_model(fmdl,inline('x<-0.4','x','y','z'));
0111    idx  = fmdl.nodes(:,1)<-0.25;
0112    fmdl.nodes(idx,1) = -0.25;
0113    fmdl.nodes(:,1) = fmdl.nodes(:,1) + 0.25;
0114    fmdl1 = crop_model(fmdl,inline('x> 1.0','x','y','z'));
0115    fmdl2 = fmdl;
0116    fmdl2.nodes(:,1) = -fmdl2.nodes(:,1);
0117    fmdl2 = crop_model(fmdl2,inline('x+y<-1.28','x','y','z'));
0118    fmdl2 = crop_model(fmdl2,inline('y<-0.95','x','y','z'));
0119    fmdl2 = crop_model(fmdl2,inline('x<-1.20','x','y','z'));
0120 show_fem(fmdl1,[0,1,1]);
0121 hold on; hh=show_fem(fmdl2,[0,1,1]); set(hh,'EdgeColor',[0,0,1]);
0122 hold off; xlim([-1.3,1.3]); axis off
0123 
0124 function [fmdl1,fmdl2]=do_unit_test_3D_mdls;
0125    imdl = mk_common_model('n3r2',[16,2]); fmdl= imdl.fwd_model;
0126    fmdl1= crop_model(fmdl,inline('x<-0.55','x','y','z'));
0127    idx  = fmdl1.nodes(:,1)<-0.25;
0128    fmdl1.nodes(idx,1) = -0.35;
0129    fmdl1.electrode([19,9])=[];
0130    fmdl1.nodes(:,1) = fmdl1.nodes(:,1) + 0.35;
0131    fmdl2 = fmdl1;
0132    fmdl2.nodes(:,1) = -fmdl2.nodes(:,1);
0133 show_fem(fmdl1,[0,1,0]);
0134 hold on; hh=show_fem(fmdl2,[0,1,0]); set(hh,'EdgeColor',[0,0,1]);
0135 hold off; xlim([-1.3,1.3]); axis off