Shape and size of anomaly does not impact reconstruction

[lrosique]

Hello,

I'm trying to use pyEIT at work for some project but I have trouble configuring it compared to EIDORS (Matlab) for an accurate size and shape of object (with dynamic jac).

When I use your example project (demo_dynamic_jac) and I change the size of the anomaly, I get this (which is perfect) : anomaly 1 = [{'x': 0.5, 'y': 0.5, 'd': 0.1, 'perm': 1000.0}] anomaly 2 = [{'x': 0.5, 'y': 0.5, 'd': 0.3, 'perm': 1000.0}] reconstruction anomaly 1 << reconstruction anomaly 2 (OK) (ie biggest anomaly = biggest activation zone)

However if I change x or y position to get closer to the center, even by a little, it's impossible to see the difference in size : anomaly 3 = [{'x': 0.4, 'y': 0.4, 'd': 0.1, 'perm': 1000.0}] anomaly 4 = [{'x': 0.4, 'y': 0.4, 'd': 0.3, 'perm': 1000.0}] reconstruction anomaly 3 == reconstruction anomaly 4 (BUG)

anomaly 5 = [{'x': 0, 'y': 0, 'd': 0.1, 'perm': 1000.0}] anomaly 6 = [{'x': 0, 'y': 0, 'd': 0.6, 'perm': 1000.0}] reconstruction anomaly 5 == reconstruction anomaly 6 (BUG)

Thank you in advance for your explanations ! and good job anyway for pyEIT, it's awesome 👍

[liubenyuan]

Hi,

if you are using eit_dynamic_jac.py, please modify the following lines

eit = jac.JAC(mesh_obj, el_pos, ex_mat=ex_mat, step=step,
              perm=1., parser='std', jac_normalized=True)
eit.setup(p=0.5, lamb=0.0001, method='kotre')
ds = eit.solve(f1.v, f0.v, normalize=True)

1, turn on jac_normalized=True 2, EIT has a low sensitivity around the center than on the boundary. This sensitivity can be weighted in the reconstruction algorithm (the p parameter in eit.setup). 3, The dynamic EIT, or time difference EIT, is a one step approximation of the static imaging. In pyeit, method='kotre' and method='lm' use the dumped least square algorithm, you can make this 1 step approximation more accurate by decreasing the penalty lambda (the lamb parameter). 4, So, you may try lamb=0.0001 for your problem. Though, smaller lamb tends to produce more ring effects around the reconstructed images. You may filter them using by spatial averging or other spatial filtering algorithm.

[lrosique]

Thank you, I've tested your jac parameters and it's better in the example but not with my dataset (even if it's slightly better). Is it possible to use DC with pyEIT ? How can we configure the input courant and voltage ? (in the code it's 1 and -1 at the injection electrodes, and changing these values yield strange results). As for using BP algorithm, i get a more accurate area but the shape is not good at all.

[liubenyuan]

1, and -1, or a and -a are current injection using simplified electrode model. What do you mean by DC (dirichlet condition?) EIDORS has many (good) reconstruction algorithms. Could you post an example data (just the current and the reference data frame, 2 frames of data in .csv format, is ok) and the expected reconstructed images using EIDORS. I can help you check whether it is the problem of pyeit.

[liubenyuan]

@lrosique Hi, your code are well writen but I have limited connection to google services.

Here are some results. Full results and the codes are emailed to you.

I add a log scale post process to pyeit.eit.jac, hope this is ok for your shape reconstruction. You should grab the latest pyeit, or manually add log_scale in your reconstructions.

The key to your reconstruction are the following codes:

# reconstruction
eit = jac.JAC(mesh_obj, elec_pos, ex_mat=ex_mat, step=step,
              perm=1., parser='std', jac_normalized=True)
eit.setup(p=0.6, lamb=0.01, method='kotre')
ds = eit.solve(f1.v, f0.v, normalize=True, log_scale=True)

1, larger p, smaller lamb may produce much better shape resolution on the center using circle anamalies, but smaller lamb may result in rings, and artefacts (spurious targets). You should find a balance tunning p, lamb 2, jac and greit (which do not use real data) are alike. This log scale post process can be applied on any reconstruction algorithms. We have used in real life applications, and it is good in finding ROI. 3, I think additional normalize on measurements row vector is not required.

[lrosique]

by DC i meant "direct current" (in EIDORS it uses Alternating Current) sorry for the abreviation ^^ for these screenshots did you use the data provided or the mesh permitivity expected in a forward to generate new data ?

[liubenyuan]

by DC i meant "direct current" (in EIDORS it uses Alternating Current) sorry for the abreviation ^^ for these screenshots did you use the data provided or the mesh permitivity expected in a forward to generate new data ?

I am using the data generated using the forward model. I tried your data, it performs bad in shape reconstruction. Is the data you provided the real measurements in physical phantoms? What shape you are using and what's the stimulation and measurement protocol?

[LU-Zirui]

Thanks for your brilliant work on pyEIT! I want to know the definition for _natural_boudary (what 1 and -1 stand for? 1A or 1V). I am using a DC power supply, so should it be a constant current or a constant voltage? and what is the unit for the conductivity in the figure?

Excuse my ignorance and looking forward to your reply.

[liubenyuan]

In Natural boundary, 1 or -1 denote 1 or -1 unit current inject on nodes. It is an ideal assumption. If you are using voltage excitation, you should modify the BC accordingly in fem.py. This unit current assumption may yield invalid reconstructed conductivities due to the mesh size.

In practice, current should flow into the domain on the edge (2D) or the face (3D) of the mesh.

[LU-Zirui]

Thank you so much for your kind reply. I'm trying to detect the conductivity distribution within a conductive-plastic circle sheet with 8 electrodes. I used 10mA constant direct current and measured the voltage drop between the electrodes in adjacent mode. I modified the eit_static_jac.py file with my data and changed _natual_boundary to +0.01/-0.01A. However, the simulation result does not converge. Is there anything wrong with my measuring method?

fv = np.array([0.1887, 0.0921, 0.0565, 0.065, 0.1809,
               0.1602, 0.0666, 0.0598,0.1022, 0.1538,
               0.1526, 0.0867, 0.0996, 0.0932, 0.1436,
               0.1547,0.106, 0.0703, 0.0672, 0.1442,
               0.1742, 0.076, 0.0527, 0.0763,0.1893,
               0.1381, 0.0621, 0.0667, 0.1047, 0.1907,
               0.1745, 0.1185,0.1263, 0.1361, 0.183,
               0.1906, 0.124, 0.0955, 0.0823, 0.1522]) # in volt

eit = jac.JAC(mesh_obj, el_pos, ex_mat, step, perm=1.0, parser="std")
eit.setup(p=0.25, lamb=1.0, method="lm")
ds = eit.gn(fv, lamb_decay=0.1, lamb_min=1e-5, maxiter=20, verbose=True)

[liubenyuan]

You can open a new issue on this.

You should pass parser=fmmu if you are using rotate-meas protocol. i.e., If your electrodes are numbered as [1..8], if now (the stimulation pattern) the the source and sink of current are 4-5, then the 5 voltgaes differences are V67, V78, V81, V12, V23.

The following are the codes, you may try using interpolation for a smoothed display

# coding: utf-8
""" demo on static solving using JAC (experimental) """
# Copyright (c) Benyuan Liu. All Rights Reserved.
# Distributed under the (new) BSD License. See LICENSE.txt for more info.
from __future__ import division, absolute_import, print_function

import numpy as np
import matplotlib.pyplot as plt

# pyEIT 2D algorithms modules
from pyeit.mesh import create
from pyeit.eit.utils import eit_scan_lines
import pyeit.eit.jac as jac

""" 1. setup """
n_el = 8
mesh_obj, el_pos = create(n_el, h0=0.12)
el_dist, step = 1, 1
ex_mat = eit_scan_lines(n_el, el_dist)
pts = mesh_obj["node"]
tri = mesh_obj["element"]

fv = np.array([0.1887, 0.0921, 0.0565, 0.065, 0.1809,
               0.1602, 0.0666, 0.0598,0.1022, 0.1538,
               0.1526, 0.0867, 0.0996, 0.0932, 0.1436,
               0.1547,0.106, 0.0703, 0.0672, 0.1442,
               0.1742, 0.076, 0.0527, 0.0763,0.1893,
               0.1381, 0.0621, 0.0667, 0.1047, 0.1907,
               0.1745, 0.1185,0.1263, 0.1361, 0.183,
               0.1906, 0.124, 0.0955, 0.0823, 0.1522]) # in volt

""" 3. solve_eit using gaussian-newton (with regularization) """
# number of stimulation lines/patterns
eit = jac.JAC(mesh_obj, el_pos, ex_mat, step, perm=1.0, parser="fmmu")
eit.setup(p=0.25, lamb=1.0, method="lm")
# lamb = lamb * lamb_decay
ds = eit.gn(fv, lamb_decay=0.1, lamb_min=1e-5, maxiter=20, verbose=True)

# plot
fig, ax = plt.subplots(figsize=(6, 4))
im = ax.tripcolor(
    pts[:, 0],
    pts[:, 1],
    tri,
    np.real(ds),
    shading="flat",
    alpha=1.0,
    cmap=plt.cm.viridis,
)
fig.colorbar(im)
ax.axis("equal")
ax.set_title("Conductivities Reconstructed")
# fig.savefig('../figs/demo_static.png', dpi=96)
plt.show()

[liubenyuan]

And you can implement your protocol and boundary conditions in fem.py