Gegularization decouplings 1ΒΆ

import crtomo
import numpy as np
from shapely.geometry import Polygon

Define the electrode positions and boundaries See Creating meshes with the mesh_gen object for more information.

electrodes = np.array((
    (0, 0),
    (1, 0),
    (2, 0),
    (3, 0),
    (4, 0),
    (5, 0),
    (6, 0),
))
boundaries = np.array((
    (-1, 0, 12),
    (0, 0, 12),
    (1, 0, 12),
    (2, 0, 12),
    (3, 0, 12),
    (4, 0, 12),
    (5, 0, 12),
    (6, 0, 12),
    (7, 0, 12),
    (7, -2, 11),
    (-1, -2, 11),
))

Additionally, we define one extra line (x0, y0, x1, y1):

extra_lines = [
    [-0, -0.20, 3, -0.45],
]

# Throws error: TODO
# poly1 = Polygon((
#    (0, 0.0),
#    (5, -0.75),
#    (5, -3),
#    (0, -3),
#    #(-2, 0),
# ))
poly1 = Polygon((
    (0, -0.20),
    (3, -0.45),
    (3, -2),
    (-1, -2),
))

"""
mesh_d = crtomo.mesh_gen.gen_mesh(
    boundaries,
    electrodes,
    # vertical line shoudl not start on horizontal line
    # additional_lines=extra_lines,
    # polygons=[poly1],
    char_lengths=[0.25, 0.5, 0.5, 0.5],
)
"""

Generate the mesh

mesh_d, _ = crtomo.mesh_gen.gen_mesh_with_polygons(
    boundaries,
    electrodes,
    # vertical line should not start on horizontal line
    # additional_lines=extra_lines,
    polygons=[poly1],
    char_lengths=[0.15, 0.25, 0.25, 0.25],

)

print(mesh_d)
mesh_d.plot_grid()
plot 12 decouplings 1
Calling the mesh generation routines...
   done
This grid was sorted using CutMcK. The nodes were resorted!
Triangular grid found
CRMod/CTRomo grid instance
number of elements: 898
number of nodes: 501
number of electrodes: 7
mesh dimensions:
X: -1.0 7.0
Z: -2.0 0.0
decouplings = crtomo.mesh_decouplings.get_decouplings(
    mesh_d, extra_lines
)
tdm = crtomo.tdMan(grid=mesh_d)
ab = tdm.configs.gen_all_current_dipoles()
tdm.configs.gen_all_voltages_for_injections(ab)
pid_mag, _ = tdm.add_homogeneous_model(300)
tdm.parman.modify_polygon(
    pid_mag,
    poly1,
    30,
)
tdm.plot_forward_models()
plot 12 decouplings 1
rmag = tdm.measurements(silent=True)[:, 0]
np.random.seed(2048)

# absolute component in [Ohm ]
noise_level_rmag_absolute = 0.01
# relative component [0, 1]
noise_level_rmag_relative = 0.05

noise_rmag = np.random.normal(
    loc=0,
    scale=rmag * noise_level_rmag_relative + noise_level_rmag_absolute
)

rmag_with_noise = rmag + noise_rmag

tdm.register_measurements(rmag_with_noise)

tdm.remove_negative_resistance_measurements()
Deleting configurations:
[]
tdm.crtomo_cfg['robust_inv'] = 'F'
tdm.add_to_decouplings(decouplings)

tdm.invert(cores=4, catch_output=True)
Attempting inversion in directory: /tmp/tmpe3t1kaza
Using binary: /usr/bin/CRTomo_dev
Calling CRTomo
Inversion attempt finished
Attempting to import the results
Reading inversion results
is robust False
Info: res_m.diag not found: /tmp/tmpe3t1kaza/inv/res_m.diag
/home/runner/work/crtomo_tools/crtomo_tools/examples/02_simple_inversion
Info: ata.diag not found: /tmp/tmpe3t1kaza/inv/ata.diag
/home/runner/work/crtomo_tools/crtomo_tools/examples/02_simple_inversion
Statistics of last iteration:
iteration                4
main_iteration           4
it_type              DC/IP
type                  main
dataRMS               2.17
magRMS                2.17
phaRMS                 NaN
lambda               91.99
roughness            5.433
cgsteps                8.0
nrdata               209.0
steplength           0.001
stepsize          0.000582
l1ratio                NaN
Name: 19, dtype: object
_ = tdm.plot_inversion_result_rmag()
plot 12 decouplings 1

Gallery generated by Sphinx-Gallery