EIDORS: Electrical Impedance Tomography and Diffuse Optical Tomography Reconstruction Software

Hosted by

Iterative Gauss Newton reconstruction in 3D

% Create simulation data $Id: basic_iterative01.m 1535 2008-07-26 15:36:27Z aadler $

% 3D Model
imdl_3d= mk_common_model('n3r2',16);
show_fem(imdl_3d.fwd_model);

sim_img= eidors_obj('image', 'stimulation image');
sim_img.fwd_model= imdl_3d.fwd_model;

% set homogeneous conductivity and simulate
sim_img.elem_data= ones( size(sim_img.fwd_model.elems,1) ,1);
homg_data=fwd_solve( sim_img );

% set inhomogeneous conductivity and simulate
sim_img.elem_data([390,391,393,396,402,478,479,480,484,486, ...
                   664,665,666,667,668,670,671,672,676,677, ...
                   678,755,760,761])= 1.15;
sim_img.elem_data([318,319,321,324,330,439,440,441,445,447, ...
                   592,593,594,595,596,598,599,600,604,605, ...
                   606,716,721,722])= 0.8;
inh_data=fwd_solve( sim_img );


slice_posn = [inf,inf,2.2,1,1; ...
              inf,inf,1.5,1,2;
              inf,inf,0.8,1,3];
show_slices(sim_img,slice_posn);

print -r75 -dpng basic_iterative01a.png

Figure: Three slices of a simple 3D shape to image (from top to bottom, at z=2.2, z=1.5, z=0.8)

Reconstruction with different iterations

% Reconstruct images $Id: basic_iterative02.m 1535 2008-07-26 15:36:27Z aadler $

% Set reconstruction parameters
imdl_3d.solve= @np_inv_solve;
imdl_3d.RtR_prior= @np_calc_image_prior;
imdl_3d.np_calc_image_prior.parameters= [3 1];
imdl_3d.hyperparameter.value= .01;

imdl_3d.fwd_model.solve=      @np_fwd_solve;
imdl_3d.fwd_model.jacobian=   @np_calc_jacobian;
imdl_3d.fwd_model.system_mat= @np_calc_system_mat;

% Number of iterations and tolerance (defaults)
imdl_3d.parameters.max_iterations = 1;
imdl_3d.parameters.term_tolerance = 1e-3;

%Add 30dB SNR noise to data
noise_level= std(inh_data.meas - homg_data.meas)/10^(30/20);
noise_level=0;
inh_data.meas = inh_data.meas + noise_level* ...
                randn(size(inh_data.meas));

% Reconstruct Images: 1 Iteration
subplot(131)
imdl_3d.parameters.max_iterations = 1;
rec_img= inv_solve(imdl_3d, homg_data, inh_data);
show_slices(rec_img,slice_posn);


% Reconstruct Images: 2 Iterations
subplot(132)
imdl_3d.parameters.max_iterations = 2;
rec_img= inv_solve(imdl_3d, homg_data, inh_data);
show_slices(rec_img,slice_posn);

% Reconstruct Images: 6 Iterations
subplot(133)
imdl_3d.parameters.max_iterations = 6;
rec_img= inv_solve(imdl_3d, homg_data, inh_data);
show_slices(rec_img,slice_posn);


print -r75 -dpng basic_iterative02a.png

Figure: Images from GN reconstructions. From left to right: 1 iteration, 2 iterations, 6 iterations. Little difference is seen in this case, mostly because this is a difference imaging problem with small contrasts.

Error norm

>>eidors_msg('log_level',3);
EIDORS:[ iter=2, norm(err)= 0.005008 ]
EIDORS:[ iter=3, norm(err)= 0.002842 ]
EIDORS:[ iter=4, norm(err)= 0.002186 ]
EIDORS:[ iter=5, norm(err)= 0.001872 ]
EIDORS:[ iter=6, norm(err)= 0.001687 ]