SIMULATE_MOVEMENT simulate small conductivity perturbations
[vh,vi,xyzr, c2f]= simulate_movement( img, xyzr );
Simulate small conductivity perturbations (using the
Jacobian calculation to be fast).
INPUT:
2D Models: specify xyzr = matrix of 3xN.
Each perturbation (i) is at (x,y) = xyzr(1:2,i)
with radius xyzr(3,i).
3D Models: specify xyzr = matrix of 4xN.
Each perturbation (i) is at (x,y,z) = xyzr(1:3,i)
with radius xyzr(4,i).
xyzr = scalar =N - single spiral of N in medium centre
img = eidors image object (with img.fwd_model FEM model).
value = the conductivity perturbation of the target (either scalar, or
vector)
OUTPUT:
vh - homogeneous measurements M x 1
vi - target simulations M x n_points
xyzr - x y and radius of each point 3 x n_points
c2f - image representation of the simulations
Example: Simulate small 3D object in centre:
img = mk_image( mk_common_model('b3cr',16) ); % 2D Image
[vh,vi] = simulate_movement( img, [0;0;0;0.05]);
Example: Simulate small 2D object in at x near side:
img = mk_image( mk_common_model('d2d3c',16) ); % 2D Image
[vh,vi] = simulate_movement( img, [0.9;0;0.05]);
show_fem(inv_solve(mk_common_model('c2c2',16),vh,vi)) % Reconst (new mesh)0001 function [vh,vi,xyzr,c2f]= simulate_movement( img, xyzr, value ); 0002 % SIMULATE_MOVEMENT simulate small conductivity perturbations 0003 % [vh,vi,xyzr, c2f]= simulate_movement( img, xyzr ); 0004 % 0005 % Simulate small conductivity perturbations (using the 0006 % Jacobian calculation to be fast). 0007 % 0008 % INPUT: 0009 % 2D Models: specify xyzr = matrix of 3xN. 0010 % Each perturbation (i) is at (x,y) = xyzr(1:2,i) 0011 % with radius xyzr(3,i). 0012 % 0013 % 3D Models: specify xyzr = matrix of 4xN. 0014 % Each perturbation (i) is at (x,y,z) = xyzr(1:3,i) 0015 % with radius xyzr(4,i). 0016 % 0017 % xyzr = scalar =N - single spiral of N in medium centre 0018 % 0019 % img = eidors image object (with img.fwd_model FEM model). 0020 % 0021 % value = the conductivity perturbation of the target (either scalar, or 0022 % vector) 0023 % 0024 % OUTPUT: 0025 % vh - homogeneous measurements M x 1 0026 % vi - target simulations M x n_points 0027 % xyzr - x y and radius of each point 3 x n_points 0028 % c2f - image representation of the simulations 0029 % 0030 % Example: Simulate small 3D object in centre: 0031 % img = mk_image( mk_common_model('b3cr',16) ); % 2D Image 0032 % [vh,vi] = simulate_movement( img, [0;0;0;0.05]); 0033 % Example: Simulate small 2D object in at x near side: 0034 % img = mk_image( mk_common_model('d2d3c',16) ); % 2D Image 0035 % [vh,vi] = simulate_movement( img, [0.9;0;0.05]); 0036 % show_fem(inv_solve(mk_common_model('c2c2',16),vh,vi)) % Reconst (new mesh) 0037 % 0038 0039 % (C) 2009 Andy Adler. Licensed under GPL v2 or v3 0040 % $Id: simulate_movement.m 4590 2014-05-02 20:46:17Z bgrychtol $ 0041 0042 if size(xyzr) == [1,1] 0043 path = linspace(0,1,xyzr); phi = 2*pi*path; 0044 meanodes= mean( img.fwd_model.nodes ); 0045 lennodes= size( img.fwd_model.nodes,1); 0046 img.fwd_model.nodes = img.fwd_model.nodes - ones(lennodes,1)*meanodes; 0047 maxnodes= max(max(abs( img.fwd_model.nodes(:,1:2) ))); 0048 img.fwd_model.nodes = img.fwd_model.nodes / maxnodes; 0049 0050 xyzr = [0.9*path.*sin(phi);0.9*path.*cos(phi);0*path; 0.05 + 0*path]; 0051 end 0052 0053 Nt = size(xyzr,2); 0054 [c2f failed] = mk_c2f_circ_mapping( img.fwd_model, xyzr); 0055 xyzr(:,failed) = []; 0056 c2f(:,failed) = []; 0057 Nt = Nt - nnz(failed); 0058 img.fwd_model.coarse2fine = c2f; 0059 0060 % We don't want a normalized jacobian here 0061 img.fwd_model = mdl_normalize(img.fwd_model,0); 0062 0063 J= calc_jacobian( img ); 0064 J= move_jacobian_postprocess( J, img, Nt); 0065 0066 vh= fwd_solve(img); 0067 vh=vh.meas; 0068 0069 vi= vh*ones(1,Nt) + J; 0070 0071 % Would this be the slow approach?: 0072 % vi = vh*zeros(1,Nt); 0073 % for i = 1: Nt 0074 % img.elem_data = 1 - c2f(:,i); 0075 % jnk = fwd_solve(img); 0076 % vi(:,i) = jnk.meas; 0077 % end 0078 0079 function J= move_jacobian_postprocess( J, img, Nt) 0080 if size(J,2) == Nt; % No problem 0081 return ; 0082 % Check if movement jacobian introduced electrode movements (elecs * dims) 0083 elseif size(J,2) == Nt + ... 0084 length(img.fwd_model.electrode) * size(img.fwd_model.nodes,2) 0085 J = J(:,1:Nt); 0086 else 0087 error('Jacobian calculator is not doing the coarse2fine mapping. This is a bug.'); 0088 end 0089 0090 0091 0092 0093 0094 % QUESTON: the accuracy of J will depend on how well we interpolate the 0095 % mesh. How important is this?