EIDORS: Electrical Impedance Tomography and Diffuse Optical Tomography Reconstruction Software

Hosted by

3D GREIT Reconstruction

Forward Model

%% Forward model
fmdl= ng_mk_cyl_models([3,2,.4],[16,1,2],[.1,0,.025]);
fmdl.stimulation = mk_stim_patterns(16,2,'{ad}','{ad}');

Simulation test data

%% Simulation
% homogeneous
imgh = mk_image(fmdl,1);
vh = fwd_solve(imgh);
% spherical target
select_fun = inline('(x-.5).^2+(y-.5).^2+(z-1.25).^2<=0.1^2','x','y','z');
% inhomogeneous
imgi = imgh;
imgi.elem_data = imgh.elem_data + elem_select(fmdl,select_fun);
vi = fwd_solve(imgi);
show_fem(imgi)
print_convert vox_GREIT_sim_02a.png

clf
show_3d_slices(imgi,1.25,.5,.5);% cuts through the target center
print_convert vox_GREIT_sim_02b.png

clf
levels(:,3) = 2.25:-.5:.25; % cuts through the target center
levels(:,1:2) = Inf;
show_slices(imgi,levels)
print_convert vox_GREIT_sim_02c.png

Figure: Three different views through a simulation target

Create 3D GREIT Model

%% Set up dual-model
imdl = select_imdl(fmdl,{'Basic GN dif'});
opt = struct();
opt.imgsz = [32 32 3*32];
opt.cube_voxels= true;
opt.prune_model= false; % % GREIT doesn't need a c2f and this is slow
imdl = mk_voxel_volume(imdl,opt);

% prune manually - remove elements with all nodes outside radius
rmdl = imdl.rec_model;
D = sqrt(sum(rmdl.nodes(:,1:2).^2,2)) <= 2;
E = all(D(reshape(rmdl.elems',4*6,[])'),2);
E2= logical(kron(E,true(6,1)));
rmdl.elems = rmdl.elems(E2,:);
rmdl.coarse2fine = rmdl.coarse2fine(E2,E);
rmdl = rmfield(rmdl, 'boundary'); 
rmdl.inside = E2;
imdl.rec_model = rmdl;

%% Set up GREIT
% specify the target distribution
[x, y, z] = ndgrid(linspace(-2,2,20),linspace(-2,2,20),.5:.25:2.5);
idx = ( x(:).^2 + y(:).^2 ) < 1.8^2;
gopt.distr = [x(idx) y(idx) z(idx)]';
% specify noise figure
gopt.noise_figure = 1.0; % if commented out, the given weight is used
% train GREIT
imdl = mk_GREIT_model(imdl,...
                      0.15,... % target radius
                      1,   ... % weight 20.5612 for NF = 1
                      gopt);   % options

%% Reconstruct
rimg = inv_solve(imdl,vh,vi);

Reconstruct and show images

show_fem(rimg)
print_convert vox_GREIT_sim_04a.png

clf
show_3d_slices(rimg,1.25,.5,.5);% cuts through the target center
print_convert vox_GREIT_sim_04b.png

clf
levels(:,3) = 2.5:-.5:.5; % cuts through the middle of the voxels
levels(:,1:2) = Inf;
show_slices(rimg,levels)
print_convert vox_GREIT_sim_04c.png

Figure: Three different views through a reconstructed object