1.1.2. tomotok.core.inversions package

Contains inversion algorithms. It is independent on the rest of the package.

Currently implemented

• Minimum Fisher Regularisation (MFR)
• Linear Algebraic Methods (LAME)
• Biorthogonal Basis Decomposition (BOB)
• Fixed parameter MFR

1.1.2.1. Submodules

1.1.2.2. tomotok.core.inversions.bob module

Contains inversion class for Biorthogonal Basis Decomposition Algorithm proposed by J. Cavlier. It is a simplified form of wavelet-vaguelette decomposition algorithm by R. Nguyen van Yen

[BOB1]

Jordan Cavalier et al., Nucl. Fusion 59 (2019): 056025

[BOB2]

Nguyen van Yen et al., Nucl. Fusion 52 (2011): 013005

class tomotok.core.inversions.bob.Bob(dec_mat=None, basis=None)

Bases: object

BiOrthogonal Basis decomposition

Attributes

basis	(scipy.sparse.dia_matrix) \(\mathbf{b}_i\) basis vectors of reconstruction plane
dec_mat	(scipy.sparse.csr_matrix) \(\hat{\mathbf{e}}_i\) decomposed matrix used to transform image into reconstruction plane
dec_mat_normed	(scipy.sparse.csr_matrix) normalised decomposed matrix
norms	(numpy.ndarray) node norms used in thresholding

Methods

decompose(gmat, basis, reg_factor=0, solver_kw: dict = None)

Decomposes the geometry matrix using basis vectors

Parameters:

gmat : scipy.sparse.csr_matrix

geometry/contribution matrix

basis : sparse matrix

matrix with basis vectors

reg_factor : float, optional

regularisation factor passed to cholesky decomposition determines weight of regularisation by identity matrix relatively to arbitrary matrix maximum value

solver_kw : dict

keyword parameters passed to the solver function

load_decomposition(floc)

Loads decomposed matrix and basis from an HDF file.

Norms are loaded only if available in file.

Parameters:

floc : str or pathlib.Path

location of hdf file with saved decomposition

normalise(precision=1e-06)

Computes normalised decomposition matrix.

Parameters:

precision : float, optional

neglects decomposition matrix rows with lower norm, by default 1e-6

save_decomposition(floc, description=”)

Saves decomposition matrix and basis to hdf file. Norms are also included if calculated.

Parameters:

floc : str or pathlib.Path

file location with name

description : str, optional

short user description for file identification

thresholding(image, c: int, precision: float = 1e-06, conv: float = 1e-09)

Applies thresholding method to provided image.

Parameters:

image : numpy.ndarray

flattened image with shape (#pixels, 1)

c : int

thresholding sensitivity constant

precision : float, optional

normalisation precision, by default 1e-6

conv : float, optional

thresholding convergence limit

Returns:

numpy.ndarray

Raises:

RuntimeError

If thresholding is called before decomposition of geometry matrix

class tomotok.core.inversions.bob.CholmodBob(dec_mat=None, basis=None)

Bases: tomotok.core.inversions.bob.Bob

Decomposition optimized for sparse matrices using Cholesky decomposition

Uses sksparse.cholmod.cholesky to solve the decomposition Requires positive definite symmetrized geometry matrix in reconstruction plane basis.

Methods

decompose(gmat, basis, reg_factor=0.001, solver_kw=None)

Decomposes geometry matrix using Cholesky decomposition and projects images

Parameters:

gmat : scipy.sparse.csr_matrix

geometry/contribution matrix

basis : sparse matrix

matrix with basis vectors

reg_factor : float, optional

regularisation factor passed to cholesky decomposition determines weight of regularisation by identity matrix relatively to arbitrary matrix maximum value

solver_kw : dict

keyword parameters passed to the solver function

class tomotok.core.inversions.bob.SimpleBob(dec_mat=None, basis=None)

Bases: tomotok.core.inversions.bob.Bob

Automatically creates simple one node basis as an sparse identity matrix

Deprecated since version 1.1.

Methods

decompose(gmat, basis=None)

class tomotok.core.inversions.bob.SparseBob(dec_mat=None, basis=None)

Bases: tomotok.core.inversions.bob.Bob

Biorthogonal Basis Decomposition optimized for sparse matrices using inverse matrix calculation.

Methods

decompose(gmat, basis, reg_factor=0, solver_kw=None)

1.1.2.3. tomotok.core.inversions.fixed module

Fixed parameter Tikhonov scheme with MFR like regularisation matrix

class tomotok.core.inversions.fixed.CholmodFixt

Bases: tomotok.core.inversions.fixed.Fixt

Cholmod version of fixed parameter MFR Requires scikit sparse to be installed in order to initialize properly.

Uses sksparse.cholmod.cholesky to solve the inversion problem

Methods

solve(a, b)

Finds solution of \(\mathbf{Ax}=\mathbf{b}\) using sksparse.cholmod.cholesky

Parameters:

a : scipy.sparse.csr_matrix

square and positive definite matrix

b : array_like

right hand side vector

class tomotok.core.inversions.fixed.Fixt

Bases: object

Inverses provided data using fixed parameter value Minimum Fisher Regularization scheme.

Attributes

last_chi	(float) person test result of the last computed inversion
logalpha	(float) \(\log(\alpha)\) logarithm of regularisation parameter
gmat	(numpy.ndarray) geometry matrix with shape (#channels, #nodes)
signal	(numpy.ndarray)
gdg	(numpy.ndarray) \(\mathbf{T}^T \cdot \mathbf{T}\) symmetrised geometry matrix, shape (#nodes, #nodes)
gdsig	(numpy.ndarray) right side of regularisation scheme \(\mathbf{T}^T \cdot \mathbf{f}\)

Methods

invert(signals, gmat, derivatives, parameters, w_factor=None, mfi_num=3, w_max=1, aniso=0)

Inverses normalised signals using fixed value of regularisation parameter.

See regularisation_matrix documentation for more information about derivative matrix formats.

Parameters:

signals : numpy.ndarray

error normalised signals on detector channels

gmat : scipy.sparse.csr_matrix

geometry matrix normalised by estimated errors

derivatives : list

list of tuples containing pairs of sparse derivatives matrices

parameters : float or list of float

regularisation parameter value(s), list length must match mfi_num

w_factor : numpy.ndarray, optional

weight matrix multipliers with shape (#nodes, )

w_max : float, optional

value used in weigh matrix for zero or negative nodes

aniso : float, optional

Determines anisotropy of derivative matrices

Returns:

g : numpy.ndarray

vector with nodes emissivities, shape (#nodes, )

stats : dict

inversion statistics

chi : float

Pearson test value of final result

niter : list of int

numbers of iterations taken in each MFI loop

logalpha : float

logarithm of final regularisation parameter

See also

regularisation_matrix

regularisation_matrix(derivatives, w, aniso=0)

Computes nonlinear regularisation matrix from provided derivative matrices and node weight factors determined by emissivity from previous iteration of the inversion loop.

Anisotropic coefficients are computed using sigmoid function.

Multiple pairs of derivatives matrices computed by different numerical scheme can be used. Each pair should contain derivatives in each locally orthogonal direction.

Parameters:

w : scipy.sparse.dia.dia_matrix

(#nodes, #nodes) diagonal matrix with pixel weight factors

derivatives : list of scipy.sparse.csc_matrix pairs

contains sparse derivative matrices, each pair contains derivatives in both locally orthogonal coordinates

aniso : float

anistropic factor, positive values make derivatives along first coordinate more significant

Returns:

scipy.sparse.csc_matrix

smoothing_mat(w, derivatives, danis)

Deprecated since version 1.1: Use regularisation_matrix()

static solve(a, b)

Finds solution of \(\mathbf{Ax}=\mathbf{b}\) using scipy.sparse.linalg.spsolve

1.1.2.4. tomotok.core.inversions.lame module

Structure of classes is based on algorithms proposed by T. Odstrcil however without sparse optimization

T. Odstrcil et al., “Optimized tomography methods for plasma emissivity reconstruction at the ASDEX Upgrade tokamak,” Rev. Sci. Instrum., 87(12), 123505.

class tomotok.core.inversions.lame.Algebraic

Bases: object

A base class for algebraic inversion methods using linear regularisation.

Attributes

u	(numpy.ndarray) decomposition matrix with shape (#channels, #channels)
s	(numpy.ndarray) diagonal from a decomposition matrix S with shape (#channels, )
v	(numpy.ndarray) decomposition matrix with shape (#nodes, #channels)

Methods

decompose(gmat, regularisation)

Prepares matrices in standard form for series expansion

Returns:u, s, v : numpy.ndarray

find_alpha(*args, **kwargs)

Finds regularisation parameter.

Returns:

float

regularisation parameter value

invert(data, gmat, regularisation, method, num=None)

Computes linear inversion using algebraic method. The inversion comprises of three stages:

• decomposition (presolving) using only geometry and derivative matrices
• searching for regularisation parameter
• solving inversion using series expansion

Parameters:

data : numpy.ndarray

gmat : numpy.ndarray

regularisation : numpy.ndarray

regularisation matrix

method : str, optional

method used for regularization parameter computation, see method find_alpha

num : int, optional

use only num most significant vectors in series expansion

Returns:

numpy.ndarray

reconstructed emissivity vector with shape (#pix,)

normalize_gmat(gmat: numpy.ndarray, errors: numpy.ndarray) → numpy.ndarray

Normalizes gmat using estimated errors

Parameters:

gmat : np.ndarray

geometry matrix

errors : np.ndarray

errors for given time slice shape (#channels,)

Returns:

np.ndarray

normalized geometry matrix

presolve(gmat, deriv)

Deprecated since version 1.1: Use decompose()

regularisation_matrix(derivatives, aniso=1)

Computes regularisation matrix from derivatives matrices.

Parameters:

derivatives : _type_

_description_

aniso : int, optional

_description_, by default 1

Returns:

_type_

_description_

series_expansion(alpha, signal, num=None)

Computes emissivity \(g\) from decomposed vectors using

\[\mathbf{g}(\alpha) = \sum_{i=1}^{m} \frac{k_{i} (\alpha)}{S_{ii}} \left( \mathbf{U}^T \cdot \mathbf{f} \cdot \tilde{\mathbf{V}} \right) {}_{*i},\]

where \(k_i(\alpha)\) are so called filtering factors computed using following formula

\[k_{i}(\alpha) = \left(1 + \frac{\alpha}{S_{ii}^2} \right)^{-1}\]

Parameters:

alpha : float

regularisation parameter

signal : numpy.ndarray

vector with channel signal

num : int

number of columns used for series expansion

Returns:

numpy.ndarray

results of inversion

class tomotok.core.inversions.lame.FastAlgebraic

Bases: tomotok.core.inversions.lame.Algebraic

A base class for linear algebraic methods using fast regularisation parameter estimation.

The regularisation parameter estimate is based on the diagonal matrix obtained by decomposition.

Methods

decompose(gmat, regularisation)

find_alpha(method: str = ‘quantile’)

Finds regularisation parameter using linear estimate based on values of diagonal.

Parameters:

method : str {‘mean’, ‘half’, ‘median’, ‘quantile’, ‘logmean}

selects method for finding regularisation parameter value, default is quantile

Returns:

float

squared value found by estimation method

class tomotok.core.inversions.lame.GevFastAlgebraic

Bases: tomotok.core.inversions.lame.FastAlgebraic

Methods

decompose(gmat, regularisation)

Decomposes geometry and regularisation matrices to form suitable for series expansion.

Uses generalised eigenvalue decomposition scheme described by L. C. Ingesson in [GEV]

Parameters:

gmat : numpy.ndarray

regularisation : numpy.ndarray

regularisation matrix with shape (#nodes, #nodes)

Returns:

u, s, v

References

[GEV]

(1, 2) L.C. Ingesson, “The Mathematics of Some Tomography Algorithms Used at JET,” JET Joint Undertaking, 2000

class tomotok.core.inversions.lame.QrFastAlgebraic

Bases: tomotok.core.inversions.lame.FastAlgebraic

Methods

decompose(gmat, regularisation)

class tomotok.core.inversions.lame.SvdFastAlgebraic

Bases: tomotok.core.inversions.lame.FastAlgebraic

Methods

decompose(gmat, regularisation)

1.1.2.5. tomotok.core.inversions.mfr module

Minimum Fisher Regularisation implemented based on articles by M. Anton and J. Mlynar

[MFR1]

Anton et al., “X-ray tomography on the TCV tokamak.”, Plasma Phys. Control. Fusion 38.11 (1996): 1849.

[MFR2]

Mlynar et al., “Current research into applications of tomography for fusion diagnostics.” J. Fusion Energy 38.3 (2019): 458-466

class tomotok.core.inversions.mfr.CholmodMfr

Bases: tomotok.core.inversions.mfr.Mfr

Uses sparse cholesky decomposition for solving the parameter optimisation task in MFI loop. Requires scikit sparse to be installed in order to initialize properly.

Uses sksparse.cholmod.cholesky to solve the regularised task in parameter optimisation.

Methods

solve(a, b)

Finds solution of \(\mathbf{Ax}=\mathbf{b}\) using sksparse.cholmod.cholesky

Parameters:

a : scipy.sparse.csr_matrix

square and positive definite matrix

b : array_like

right hand side vector

class tomotok.core.inversions.mfr.Mfr

Bases: object

Inverses provided data using Minimum Fisher Regularization scheme.

Attributes

logalpha	(float) \(\log(\alpha)\) logarithm of regularisation parameter
gmat	(numpy.ndarray) geometry matrix with shape (#channels, #nodes)
signal	(numpy.ndarray)
gdg	(numpy.ndarray) \(\mathbf{T}^T \cdot \mathbf{T}\) symmetrised geometry matrix, shape (#nodes, #nodes)
gdsig	(numpy.ndarray) right side of regularisation scheme \(\mathbf{T}^T \cdot \mathbf{f}\)

Methods

determine_regularisation(derivatives, w, aniso, bounds, iter_max, tolerance)

Uses minimize scalar function from scipy to iteratively minimize chi square (Pearson test).

invert(signals, gmat, derivatives, w_factor=None, mfi_num=3, bounds=(-15, 0), iter_max=15, w_max=1, aniso=0, tolerance=0.05, zero_negative=False)

Inverses normalised signals using mfi_num Fisher Information cycles each with iter_max steps of regularisation parameter optimisation.

See regularisation_matrix documentation for more information about derivative matrix formats.

Parameters:

signals : numpy.ndarray

error normalised signals on detector channels

gmat : scipy.sparse.csr_matrix

geometry matrix normalised by estimated errors, shape (#channels, #nodes)

derivatives : list

list of tuples containing pairs of sparse derivatives matrices

w_factor : numpy.ndarray, optional

weight matrix multipliers with shape (#nodes, )

mfi_num : float, optional

number of MFR iterations

bounds : tuple of two floats, optional

exponent values for bounds of regularisation parameter alpha

iter_max : float, optional

maximum number of root finding iterations

w_max : float, optional

value used in weigh matrix for zero or negative nodes

aniso : float, optional

Determines anisotropy of derivative matrices

tolerance : float, optional

End MFI loop if residuum reaches interval (1-tolerance; 1+tolerance), by default 0.05

zero_negative : bool

sets negative values in previous step result to zero when initializing MFI loop

Returns:

g : numpy.ndarray

vector with nodes emissivities, shape (#nodes, )

stats : dict

inversion statistics

chi : list of float

Pearson test value of each optimization result

iter_num : list of int

numbers of iterations taken in each MFI loop

logalpha : list of float

final regularisation parameter logarithm of each MFI loop

See also

regularisation_matrix

regularisation_matrix(derivatives, w, aniso=0)

Computes nonlinear regularisation matrix from provided derivative matrices and node weight factors determined by emissivity from previous iteration of the inversion loop.

Anisotropic coefficients are computed using sigmoid function.

Multiple pairs of derivatives matrices computed by different numerical scheme can be used. Each pair should contain derivatives in each locally orthogonal direction.

Parameters:

w : scipy.sparse.dia.dia_matrix

(#nodes, #nodes) diagonal matrix with pixel weight factors

derivatives : list of scipy.sparse.csc_matrix pairs

contains sparse derivative matrices, each pair contains derivatives in both locally orthogonal coordinates

aniso : float

anistropic factor, positive values make derivatives along first coordinate more significant

Returns:

scipy.sparse.csc_matrix

smoothing_mat(w, derivatives, danis)

Deprecated since version 1.1: Use regularisation_matrix()

static solve(a, b)

Finds solution of \(\mathbf{Ax}=\mathbf{b}\) using scipy.sparse.linalg.spsolve