EIDORS: Electrical Impedance Tomography and Diffuse Optical Tomography Reconstruction Software

Hosted by

Compare 2D image reconstructions

% Compare 2D algorithms
% $Id: tutorial120a.m 3273 2012-06-30 18:00:35Z aadler $

imb=  mk_common_model('c2c',16);

e= size(imb.fwd_model.elems,1);
bkgnd= 1;

% Solve Homogeneous model
img= mk_image(imb.fwd_model, bkgnd);
vh= fwd_solve( img );

% Add Two triangular elements
img.elem_data([25,37,49:50,65:66,81:83,101:103,121:124])=bkgnd * 2;
img.elem_data([95,98:100,79,80,76,63,64,60,48,45,36,33,22])=bkgnd * 2;
vi= fwd_solve( img );

% Add -12dB SNR
vi_n= vi; 
nampl= std(vi.meas - vh.meas)*10^(-18/20);
vi_n.meas = vi.meas + nampl *randn(size(vi.meas));

show_fem(img);
axis square; axis off
print_convert('tutorial120a.png', '-density 60')

Figure: Sample image to test 2D image reconstruction algorithms

% Compare 2D algorithms
% $Id: tutorial120b.m 5522 2017-06-07 12:03:37Z aadler $
clf;clear imgr imgn

% Create Inverse Model
inv2d= eidors_obj('inv_model', 'EIT inverse');
inv2d.reconst_type= 'difference';
inv2d.jacobian_bkgnd.value= 1;

% This is not an inverse crime; inv_mdl != fwd_mdl
imb=  mk_common_model('b2c',16);
inv2d.fwd_model= imb.fwd_model;

% Guass-Newton solvers
inv2d.solve=       @inv_solve_diff_GN_one_step;

% Tikhonov prior
inv2d.hyperparameter.value = .03;
inv2d.RtR_prior=   @prior_tikhonov;
imgr(1)= inv_solve( inv2d, vh, vi);
imgn(1)= inv_solve( inv2d, vh, vi_n);

% NOSER prior
inv2d.hyperparameter.value = .1;
inv2d.RtR_prior=   @prior_noser;
imgr(2)= inv_solve( inv2d, vh, vi);
imgn(2)= inv_solve( inv2d, vh, vi_n);

% Laplace image prior
inv2d.hyperparameter.value = .1;
inv2d.RtR_prior=   @prior_laplace;
imgr(3)= inv_solve( inv2d, vh, vi);
imgn(3)= inv_solve( inv2d, vh, vi_n);

% Automatic hyperparameter selection
inv2d.hyperparameter = rmfield(inv2d.hyperparameter,'value');
inv2d.hyperparameter.func = @choose_noise_figure;
inv2d.hyperparameter.noise_figure= 0.5;
inv2d.hyperparameter.tgt_elems= 1:4;
inv2d.RtR_prior=   @prior_gaussian_HPF;
inv2d.solve=       @inv_solve_diff_GN_one_step;
imgr(4)= inv_solve( inv2d, vh, vi);
imgn(4)= inv_solve( inv2d, vh, vi_n);
inv2d.hyperparameter = rmfield(inv2d.hyperparameter,'func');

% Total variation using PDIPM
inv2d.hyperparameter.value = 1e-5;
inv2d.solve=       @inv_solve_TV_pdipm;
inv2d.R_prior=     @prior_TV;
inv2d.parameters.max_iterations= 10;
inv2d.parameters.term_tolerance= 1e-3;

%Vector of structs, all structs must have exact same (a) fields (b) ordering
imgr5= inv_solve( inv2d, vh, vi);
imgr5=rmfield(imgr5,'type'); imgr5.type='image';
imgr(5)=imgr5;
imgn5= inv_solve( inv2d, vh, vi_n);
imgn5=rmfield(imgn5,'type'); imgn5.type='image';
imgn(5)=imgn5;

% Output image
imgn(1).calc_colours.npoints= 128;
imgr(1).calc_colours.npoints= 128;

show_slices(imgr, [inf,inf,0,1,1]);
print_convert tutorial120b.png;

show_slices(imgn, [inf,inf,0,1,1]);
print_convert tutorial120c.png;

Figure: Images reconstructed with data without noise, From Left to Right: 1) One step Gauss-Newton reconstruction (Tikhonov prior) 2) One step Gauss-Newton reconstruction (NOSER prior) 3) One step Gauss-Newton reconstruction (Laplace filter prior) 5): One step Gauss-Newton reconstruction (automatic hyperparameter selection) 5): Total Variation reconstruction

Figure: Images reconstructed with data with added 12dB SNR. From Left to Right: 1) One step Gauss-Newton reconstruction (Tikhonov prior) 2) One step Gauss-Newton reconstruction (NOSER prior) 3) One step Gauss-Newton reconstruction (Laplace filter prior) 5): One step Gauss-Newton reconstruction (automatic hyperparameter selection) 5): Total Variation reconstruction