jacobian_3d_fields

 [J] = jacobian_3d_fields(V,Ela,D,elec,vtx,simp,mat_ref,v_f,df, c2f);

calculates the Jacobian (sensitivity) matrix from V fwd

J        = The Jacobian (sensitivity) matrix with respect to conductivity

V        = forward solver voltage
D,Ela    = parameters from the fem_master_full
elec     = the electrodes matrix
vtx      = The vertices matrix
simp     = The simplices matrix
mat_ref  = The reference conductivity vector
v_f      = The measurement fields
df       = Measurements per current pattern as used in v_f 
c2f      = Coarse to fine map between fine and coarse model (optional)

0001 function [J] = jacobian_3d_fields(V,Ela,D,elec,vtx,simp,mat_ref,v_f,df, c2f);
0002 % [J] = jacobian_3d_fields(V,Ela,D,elec,vtx,simp,mat_ref,v_f,df, c2f);
0003 %
0004 %calculates the Jacobian (sensitivity) matrix from V fwd
0005 %
0006 %J        = The Jacobian (sensitivity) matrix with respect to conductivity
0007 %
0008 %V        = forward solver voltage
0009 %D,Ela    = parameters from the fem_master_full
0010 %elec     = the electrodes matrix
0011 %vtx      = The vertices matrix
0012 %simp     = The simplices matrix
0013 %mat_ref  = The reference conductivity vector
0014 %v_f      = The measurement fields
0015 %df       = Measurements per current pattern as used in v_f
0016 %c2f      = Coarse to fine map between fine and coarse model (optional)
0017 
0018 % $Id: jacobian_3d_fields.html 2819 2011-09-07 16:43:11Z aadler $
0019 
0020 [vr,dim] = size(vtx);
0021 
0022 if sum(df)~= size(v_f,2);
0023    error('Mismatched data input');
0024 end
0025 
0026 %Select the part referring to the interior nodes
0027 V = V(1:vr,:);
0028 v_f = v_f(1:vr,:);
0029 
0030 n_elem= size(simp,1);
0031 %diag_Ela = diag(Ela(1:dim:size(Ela,1),1:dim:size(Ela,2)));
0032 diag_Ela = diag(Ela(1:dim:size(Ela,1),1:dim:size(Ela,2)));
0033 diag_Ela = spdiags(diag_Ela, 0, n_elem, n_elem);
0034 
0035 if nargin>=10 % coarse2fine provided
0036    J = zeros(sum(df), size(c2f,2) );
0037    diag_Ela = diag_Ela*c2f;
0038 else
0039    J = zeros(sum(df),n_elem );
0040 end
0041 cnt = 0;
0042 
0043 for p=1:size(V,2)  % for each stimulation pattern
0044 
0045    DV =  D*V(:,p); %Gradient of the current fields
0046    df_idx= sum(df(1:p-1));
0047 
0048    for m=1:df(p)   % for each measurement in this stim
0049 
0050      Dvf = D*v_f(:,df_idx + m); %Gradient of the measurement fields
0051 
0052      Jrow_x3 = Dvf .* DV ;  
0053      Jrow_u = sum(reshape(Jrow_x3,dim,[]),1); % Works for 2D and 3D now
0054      
0055 %    Jrow = Jrow_u .* diag_Ela.';
0056      Jrow = Jrow_u * diag_Ela;
0057      
0058      cnt = cnt+1;
0059      J(cnt,:) = -Jrow;
0060      
0061   end %m
0062 end %p
0063   
0064 
0065 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0066 % This is part of the EIDORS suite.
0067 % Copyright (c) N. Polydorides and W.R.B. Lionheart 2003
0068 % Copying permitted under terms of GNU GPL
0069 % See enclosed file gpl.html for details.
0070 % EIDORS 3D version 2.0
0071 % MATLAB version 5.3 R11
0072 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0073