API reference: The icsd3d package

class icsd3d.icsd3d_class_dev.iCSD3d_Class(dirName)

Create a icsd inversion object.

Parameters:

• coord_file (str, mandatory) – coordinates of the VRTe for plotting
• wr (float, optional) – Weight regularization
• wc (float, optional) –

Methods

Invert([pareto])	Invert the voltage to current densities.
con_A_f()	Set current conservation constrainst on A (rows of ones)
con_b_f()	Set current conservation constrainst on b
con_w_f()	Set current conservation constrainst weight; default is wc=1e6
createSurvey()	Data container for survey paramaters such as geometry file
createTimeLapseSurvey(fnames)	Import multiple surveys.
createdirs()
iCSD()	solve linear system, given A matrix (VRTe, constrain, regul) and b (observations)
icsd_init()	these functions are called only once, even for pareto, as they do not depend on the regularization weight
labelsM0(method)	Parse graphical labels to plot estimate M0 model
load_coord()
mkGrid_XI_YI()	grid for interpolation
normF1()
nx_ny()	find number of nodes in each direction, has to be a regular grid
nx_ny_nz()	find number of nodes in each direction, has to be a regular grid
obs_w_f()	weight the observations, can also ignore observations by setting w = 0
parseDataReg()	Parse regularisation parameters before inversion
parseM0(method)	Parse regularisation parameters before inversion
parseModelReg()	Parse regularisation parameters before inversion
plotCSD()	Plot CSD in 2d, using matplotlib and scipy interpolation
plot_knee_icsd()	Plot CSD for the best regularisation parameter after L-curve automatic analysis using a knee-locator
prepare4iCSD()	this function is called for each weight, keep them separated for pareto
regularize_A()	create and append rows for spatial regularization to A
regularize_A_3d()	create and append rows for spacial regularization to A
regularize_A_x_y()	create and append rows for spatial regularization to A, second derivative is applied in both direction x and y math:: Dx = ?? Dy= We used a ponderate diagonal matrix with coeffcient 1,-2, 1
regularize_smallnessX0()	Create relative smallness instance
regularize_sum_AX0()	sum smallness and spatial regularisation
regularize_w()	create vector with weights, the length is determined by the number of regul rows in A
run_pareto()	run iCSD multiple times while changing the weights to explore the L-curve
run_productmoment()	Compute the product moment correlation after Binley et al. 1999
run_single()	Run a single inversion (unique regularisation weight) Equivalent to several steps:: self.prepare4iCSD() self.plotCSD() self.RMSAnalysis()
saveInvData(outputdir)	Save inverted data
showResults([ax, clim, cmap, plotElecs, sc, …])	Show inverted model.
stack_A()	Stack A (green fcts), constrainsts and regularisation
stack_b()	Stack b, constrainsts and regularisation
stack_w()	create vector with weights for observation, constrain, and regularization then use it as diagonal for the weight matrix
weight_A()	Apply the weights to A
weight_b()	Apply the weights to b

DetectKneePt
Export_sol
RMSAnalysis
ResidualAnalysis
check_nVRTe
estimateM0
load_geom
load_obs
load_sim
misfit_2_initialX0
plotCSD3d
plotCSD3d_pyvista
plotScattered3d
plotmisfitF1
regularize_A_UnstructuredMesh3d
regularize_A_x_y_z
regularize_b
reshape_A
run_misfitF1
showEstimateM0
writeFIT

Invert(pareto=False)

Invert the voltage to current densities.

con_A_f()

Set current conservation constrainst on A (rows of ones)

con_b_f()

Set current conservation constrainst on b

con_w_f()

Set current conservation constrainst weight; default is wc=1e6

createSurvey()

Data container for survey paramaters such as geometry file

createTimeLapseSurvey(fnames)

Import multiple surveys.

Parameters:fnames (list of str) – List of file to be parsed or directory where the files are.

iCSD()

solve linear system, given A matrix (VRTe, constrain, regul) and b (observations)

icsd_init()

these functions are called only once, even for pareto, as they do not depend on the regularization weight

labelsM0(method)

Parse graphical labels to plot estimate M0 model

mkGrid_XI_YI()

grid for interpolation

nx_ny()

find number of nodes in each direction, has to be a regular grid

nx_ny_nz()

find number of nodes in each direction, has to be a regular grid

obs_w_f()

weight the observations, can also ignore observations by setting w = 0

parseDataReg()

Parse regularisation parameters before inversion

parseM0(method)

Parse regularisation parameters before inversion

parseModelReg()

Parse regularisation parameters before inversion

plotCSD()

Plot CSD in 2d, using matplotlib and scipy interpolation

Parameters:self –

plot_knee_icsd()

Plot CSD for the best regularisation parameter after L-curve automatic analysis using a knee-locator

Parameters:self –

prepare4iCSD()

this function is called for each weight, keep them separated for pareto

regularize_A()

create and append rows for spatial regularization to A

regularize_A_3d()

create and append rows for spacial regularization to A

regularize_A_x_y()

create and append rows for spatial regularization to A, second derivative is applied in both direction x and y math:: Dx = ?? Dy= We used a ponderate diagonal matrix with coeffcient 1,-2, 1

regularize_smallnessX0()

Create relative smallness instance

\[X_{0} = A*lpha_{x_{0}}\]

Parameters:self –

regularize_sum_AX0()

sum smallness and spatial regularisation

\[W_{m}=lpha_{s}I+{D_{x}}^{T}D_{x} + D_{z}}^{T}D_{z}\]

Parameters:self –

regularize_w()

create vector with weights, the length is determined by the number of regul rows in A

run_pareto()

run iCSD multiple times while changing the weights to explore the L-curve

run_productmoment()

Compute the product moment correlation after Binley et al. 1999

\[r_{k}=\]

rac{sum_{i}(D_{I}-overline{D})(F_{i}(I_{k})-overline{F}(I_{k}))}{sqrt{sum_{i}(D_{I}-overline{D})^{2}}sum_{i}(F_{i}(I_{k})-overline{F}(I_{k}))^{2}}

where $D_{i}$ is the $i^{th}$ measured transfer resistance and $F_{i}(I_{k})$ is the $i^{th}$ transfer resistance computed to unit current at location k.

run_single()

Run a single inversion (unique regularisation weight) Equivalent to several steps:

self.prepare4iCSD()
self.plotCSD()
self.RMSAnalysis()

saveInvData(outputdir)

Save inverted data

Parameters:outputdir (str) – Path where the .csv files will be saved.

showResults(ax=None, clim=None, cmap='viridis_r', plotElecs=False, sc=None, retElec=None, mesh=None, gif3d=False, title=None)

Show inverted model.

Parameters:

• ax (Matplotlib.Axes, optional) – If specified, the graph will be plotted against this axis.
• clim (array, optional) – Minimum and Maximum value of the colorbar.
• cmap (str, optional) – Name of the Matplotlib colormap to use.
• plotElecs (bool, optional) – If True add to the ICSD plot measuring electrodes as points
• sc (array, optional) – Coordinates of the sources, format = x1,y1 x2,y2’ If Not None add to the ICSD plot the source A electrode
• retElec (array, optional) – Coordinates of the return electrode, format = x1,y1’) If Not None add to the ICSD plot the return B electrode
• mesh (str, optional) – Specify name of the vtk file If Not None add mesh3d.vtk to plot with the results of icsd (for 3d using pyvista)
• gif3d (bool, optional) –

  If True record a gif using orbital function of pyvista title : str, optional

  Specify inversion titlename to be add to the plot

stack_A()

Stack A (green fcts), constrainsts and regularisation

stack_b()

Stack b, constrainsts and regularisation

stack_w()

create vector with weights for observation, constrain, and regularization then use it as diagonal for the weight matrix

weight_A()

Apply the weights to A

weight_b()

Apply the weights to b

icsd3d.inversion: inversion scheme

Prior model

Created on Mon May 11 16:22:08 2020 @author: Benjamin Estimation of initial model based on the physical assumption that a single source current can describe the pattern of the masse anomaly

inversion.priorM0.normF1(A, b)

compute the norm between observation data and individual green functions

inversion.priorM0.productmoment(A, b)

Compute the product moment correlation after Binley et al. 1999

\[r_{k}=\]

rac{sum_{i}(D_{I}-overline{D})(F_{i}(I_{k})-overline{F}(I_{k}))}{sqrt{sum_{i}(D_{I}-overline{D})^{2}}sum_{i}(F_{i}(I_{k})-overline{F}(I_{k}))^{2}}

where $D_{i}$ is the $i^{th}$ measured transfer resistance and $F_{i}(I_{k})$ is the $i^{th}$ transfer resistance computed to unit current at location k.

Smoothing

Created on Mon May 11 17:29:01 2020 @author: Benjamin

inversion.smoothing.nx_ny(coord)

find number of nodes in each direction, has to be a regular grid

inversion.smoothing.nx_ny_nz(coord)

find number of nodes in each direction, has to be a regular grid

inversion.smoothing.ponderate_smallnessX0(alphaSxy, alphax0, **kwargs)

Create relative smallness instance and applied smallness coefficient (lpha_{x_{0}}) weight

\[X_{0} = A*lpha_{x_{0}}\]

Parameters:self –

inversion.smoothing.regularize_A(coord, nVRTe)

create and append rows for to A, for spatial regularization (simple model smoothing). Working only on 2d regular meshes

inversion.smoothing.regularize_A_3d(nVRTe, coord)

model smoothing consisting in creating and appending rows for spatial regularization to A

inversion.smoothing.regularize_A_UnstructuredMesh3d(coord, nVRTe, k_neighbors=4)

model smoothing consisting in creating and appending rows for spatial regularization to A. Adapted for unstructured mesh since it uses the k_neighbors method, default k=4. Also working on regular grid 2d

inversion.smoothing.regularize_A_x_y(coord, alphaSx, alphaSy)

create and append rows for spatial regularization to A, second derivative is applied in both direction x and y math:: Dx = ?? Dy= We used a ponderate diagonal matrix with coeffcient (1,-2, 1)

inversion.smoothing.regularize_A_x_y_z(coord)

Model smoothing in 3d, not tested not working

inversion.smoothing.regularize_b(reg_A)

initiate vector b with zeros, the length is determined by the number of regul rows in A

inversion.smoothing.regularize_w(reg_A, wr, x0_prior, **kwargs)

create vector with weights, the length is determined by the number of regul rows in A such as .. math :: A = (G’*Wd*G + lambda*Wm)

b = G’*Wd*d + lambda*Wm*m0;

inversion.smoothing.sum_smallness_smoothness(alphaSxy, x0_prior, **kwargs)

sum smallness and spatial regularisation

\[W_{m}=lpha_{s}I+{D_{x}}^{T}D_{x} + D_{z}}^{T}D_{z}\]

Parameters:self –

Solver

Created on Tue May 12 09:35:37 2020

@author: Benjamin

inversion.solve.con_A_f(A)

Set current conservation constrainst on A (rows of ones)

inversion.solve.con_b_f(b)

Set current conservation constrainst on b

inversion.solve.con_w_f(wc)

Set current conservation constrainst weight; default is wc=1e6

inversion.solve.iCSD(x0_ini_guess, A_w, b_w, dim, coord, path, **kwargs)

Solve linear system, given weigted A matrix (VRTe, constrain, regul) and weigted b (observations).

Parameters:

• x0_ini_guess (*) – Initial guess
• A_w (*) – Kernel of green functions
• b_w (*) – Weigted observations
• dim (*) – Survey dimension i.e 2d or 3d
• coord (*) – Coordinates of the virtual sources

Returns:

x – Solution

Return type:

1D-arrays

inversion.solve.obs_w_f(obs_err, b, errRmin, sd_rec=None)

weight the observations, can also ignore observations by setting w = 0

inversion.solve.stack_A(A, con_A, reg_A)

Stack A (green fcts), constrainsts and regularisation

inversion.solve.stack_b(b, con_b, reg_b)

Stack b, constrainsts and regularisation

inversion.solve.stack_w(obs_w, con_w, x0_prior, **kwargs)

create vector with weights for observation, constrain, and regularization then use it as diagonal for the weight matrix

inversion.solve.weight_A(x0_prior, A_s, **kwargs)

Apply the weights to A

inversion.solve.weight_b(x0_prior, b_s, **kwargs)

Apply the weights to b

icsd3d.plotters: plotters for results visualisation

Created on Mon May 11 14:44:26 2020

@author: Benjamin

plotters.mpl_plot.labels(method)

Parse graphical labels to plot

plotters.mpl_plot.plotCSD2d(coord, data_sol, b, b_w, xfun, path, pareto, retElec=None, sc=None, ax=None, **kwargs)

Plot CSD in 2d, using matplotlib and scipy interpolation

Parameters:self –

plotters.mpl_plot.plotCSD3d(wr, coord, data, path, filename, knee, KneeWr, ax=None, title=None, pltRemotes=False, **kwargs)

plot scattered 3d current sources density for a given regularisation weight wr (can be the knee location if pareto-curve mode is run)

Parameters:

• sc (sources coordinates) –
• (to add) (kwargs) –

plotters.mpl_plot.plotContour2d(coord, data_sol, physLabel, path, retElec=None, sc=None, **kwargs)

Plot contour in 2d, using matplotlib and scipy interpolation

Parameters:self –

plotters.mpl_plot.plot_knee_icsd(wr, kn)

Plot CSD for the best regularisation parameter after L-curve automatic analysis using a knee-locator

Parameters:self –

plotters.mpl_plot.showObs2d(path, **kwargs)

Plot contour in 2d, using matplotlib and scipy interpolation. Required surface and borehole electrode to make the 2d interpolation possible

Parameters:self –

Created on Mon May 11 14:44:26 2020 @author: Benjamin 3D plots using pyvista

icsd3d.gridder: Grid data

Created on Mon May 11 15:58:36 2020

@author: Benjamin

gridder.mkgrid.mkGrid_XI_YI(coord_x, coord_y, nbe=500)

Grid for interpolation

Parameters:

• coord_x (*) – The x coordinates of the grid points
• coord_y (*) – The y coordinate of the grid points

Returns:

XI, YI – Meshgrid points

Return type:

1D-arrays

icsd3d.importers: wrappers to facilitate import of common ERT data

Created on Mon May 11 15:18:31 2020

@author: Benjamin

importers.read.DataImport(SimFile=None, ObsFile=None)

Data importer for common data files (Resipy and Gimli) Import and parse observation files, simulated file and geometry file

importers.read.load_coord(path, filename, dim)

load coordinates of the virtual current sources

importers.read.load_geom(path)

load the geometry of the acquisition (*geom file custum for Mise-à-la-masse data)

importers.read.load_obs(path, filename)

load the observations file (normalised voltages)

importers.read.load_sim(path, filename)

load the simulated green functions file