Surface nuclear magnetic resonance (NMR) data inversion



MRS([name, verbose])

Magnetic resonance sounding (MRS) manager class.

MRS1dBlockQTModelling(nlay, K, zvec, t[, ...])

MRS1dBlockQTModelling - pygimli modelling class for block-mono QT inversion

MRSprofile([filename, x, dx, x0])

manager class for several MRS data along a profile for joint inversion


class pygimli.physics.sNMR.MRS(name=None, verbose=True, **kwargs)[source]#

Bases: object

Magnetic resonance sounding (MRS) manager class.

  • q (t,) –

  • error (data,) –

  • z (K,) –

  • modelU (model, modelL,) –

loadMRSI - load MRSI (MRSmatlab format) data
showCube - show any data/error/misfit as data cube (over q and t)
showDataAndError - show data and error cubes
showKernel - show Kernel matrix
createFOP - create forward operator
createInv - create pygimli Inversion instance
run - run block-mono (alternatively smooth-mono) inversion (with bootstrap)
calcMCM - compute model covariance matrix and thus uncertainties
splitModel - return thickness, water content and T2* time from vector
showResult/showResultAndFit - show inversion result (with fit)
runEA - run evolutionary algorithm (GA, PSO etc.) using inspyred
plotPopulation - plot final population of an EA run
__init__(name=None, verbose=True, **kwargs)[source]#

MRS init with optional data load from mrsi file

  • name (string) – Filename with load data and kernel (.mrsi) or just data (.mrsd)

  • verbose (bool) – be verbose

:param kwargs - see MRS.loadMRSI().:


Compute linear model covariance matrix.


Compute model bounds using covariance matrix diagonals.


Check data and retrieve data and error vector.

static createFOP(nlay, K, z, t)[source]#

Create forward operator instance.

createInv(nlay=3, lam=100.0, verbose=True, **kwargs)[source]#

Create inversion instance (and fop if necessary with nlay).


Generate (GA) model from random vector (0-1) using model bounds.

invert(nlay=3, lam=100.0, startvec=None, verbose=True, uncertainty=False, **kwargs)[source]#

Easiest variant doing all (create fop and inv) in one call.


Load data cube from single ascii file (old stuff)

loadDataNPZ(filename, **kwargs)[source]#

Load data and kernel from numpy gzip packed file.

The npz file contains the fields: q, t, D, (E), z, K


Load several standard files from dir (old Borkum stage).


Load error cube from a single ascii file (old stuff).


Load kernel matrix from mrsk or two bmat files.

loadKernelNPZ(filename, **kwargs)[source]#

Load data and kernel from numpy gzip packed file.

The npz file contains the fields: q, t, D, (E), z, K

loadMRSD(filename, usereal=False, mint=0.0, maxt=2.0)[source]#

Load mrsd (MRS data) file: not really used as in MRSD.

loadMRSI(filename, **kwargs)[source]#

Load data, error and kernel from mrsi or mrsd file

  • usereal (bool [False]) – use real parts (after data rotation) instead of amplitudes

  • mint/maxt (float [0.0/2.0]) – minimum/maximum time to restrict time series


Load inversion result from column file.


Load the kernel vertical discretisation (z) vector.


Plot EA statistics (best, worst, …) over time.

plotPopulation(maxfitness=None, fitratio=1.05, savefile=True)[source]#

Plot fittest individuals (fitness<maxfitness) as 1d models

  • maxfitness (float) – maximum fitness value (absolute) OR

  • fitratio (float [1.05]) – maximum ratio to minimum fitness


Return block model results (thk, wc and T2 vectors).

run(verbose=True, uncertainty=False, **kwargs)[source]#

Easiest variant doing all (create fop and inv) in one call.

runEA(nlay=None, eatype='GA', pop_size=100, num_gen=100, runs=1, mp_num_cpus=8, **kwargs)[source]#

Run evolutionary algorithm using the inspyred library

  • nlay (int [taken from classic fop if not given]) – number of layers

  • pop_size (int [100]) – population size

  • num_gen (int [100]) – number of generations

  • runs (int [pop_size*num_gen]) – number of independent runs (with random population)

  • eatype (string ['GA']) –

    algorithm, choose among:

    ’GA’ - Genetic Algorithm [default] ‘SA’ - Simulated Annealing ‘DEA’ - Discrete Evolutionary Algorithm ‘PSO’ - Particle Swarm Optimization ‘ACS’ - Ant Colony Strategy ‘ES’ - Evolutionary Strategy

saveFigs(basename=None, extension='pdf')[source]#

Save all figures to (pdf) files.


Save inversion result to column text file for later use.


Set parameter boundaries for inversion.

showCube(ax=None, vec=None, islog=None, clim=None, clab=None)[source]#

Plot any data (or response, error, misfit) cube nicely.

showDataAndError(figsize=(10, 8), show=False)[source]#

Show data cube along with error cube.


Show the kernel as matrix (Q over z).

showResult(figsize=(10, 8), save='', fig=None, ax=None)[source]#

Show theta(z) and T2*(z) (+uncertainties if there).

showResultAndFit(figsize=(12, 10), save='', plotmisfit=False, maxdep=0, clim=None)[source]#

Show ec(z), T2*(z), data and model response.

static simulate(model, K, z, t)[source]#

Do synthetic modelling.


Split model vector into d, theta and T2*.

class pygimli.physics.sNMR.MRS1dBlockQTModelling(nlay, K, zvec, t, verbose=False)[source]#

Bases: ModellingBaseMT__

MRS1dBlockQTModelling - pygimli modelling class for block-mono QT inversion

f=MRS1dBlockQTModelling(lay, KR, KI, zvec, t, verbose = False )

__init__(nlay, K, zvec, t, verbose=False)[source]#

Constructor with number of layers, kernel, z and t vectors.


Yield model response cube as vector.

class pygimli.physics.sNMR.MRSprofile(filename=None, x=None, dx=1, x0=0, **kwargs)[source]#

Bases: object

manager class for several MRS data along a profile for joint inversion

  • mrs (list of MRS objects (single soundings)) –

  • x (list of positions for the soundings) –

load - load mrs files from a directory
set X - set x vector
showData - show MRS data
independentBlock1dInversion - perform independent 1D block inversion
block1dInversion - 1D block inversion of all data sets together
blockLCInversion - 1D block laterally constrained inversion of all data
printFits - print total misfit (chi^2, rms) and individual values
showModel - show LCI model
__init__(filename=None, x=None, dx=1, x0=0, **kwargs)[source]#

Initialize profile object by mrs objects and optional positions.

  • filename (list of str | str) – list of files OR filenames(with *) OR directory to load

  • x (iterable) – position vector of individual soundings

  • x0 (float [0]) – starting position

  • dx (position [1]) – position increment

block1dInversion(nlay=2, lam=100.0, show=False, verbose=True, uncertainty=False)[source]#

Invert all data together by a 1D model (more general solution).

block1dInversionOld(nlay=2, startModel=None, verbose=True, uncertainty=False, **kwargs)[source]#

Invert all data together by one 1D model (variant 1 - all equal).

blockLCInversion(nlay=2, startModel=None, **kwargs)[source]#

Laterally constrained (piece-wise 1D) block inversion.

independentBlock1dInversion(nlay=2, lam=100, startModel=None)[source]#

Independent inversion of all soundings.

  • nlay (int [2]) – number of layers

  • lam (float) – regularization parameter

  • startModel (array/vector) – starting model (see parameters)

load(filenames, **kwargs)[source]#

load mrs files in a list of (single) MRS handlers filename can be a list of mrsi files or a directory to search Additional parameters: usereal, mint, maxt (see MRS.load)


load one kernel file for all soundings


Show single fits and total fit.

saveFigs(basename='out', extension='pdf')[source]#

Save all figures to (pdf) files.

setX(x=None, x0=0, dx=1)[source]#

define positions for soundings and sort accordingly


Show 1D model (e.g. of joint block inversion).

showData(figsize=(15, 10), nc=0, nr=0, clim=None)[source]#

show all data cubes in subplots


Show chi-square and rms fits of individual soundings.


show initial values of whole profile

showModel(showFit=0, cmap='Spectral', figsize=(13, 12), wlim=(0, 0.5), tlim=(0.05, 0.5))[source]#

Show 2d model as stitched 1d models along with fit.

showT2(tlim=(0.05, 0.5), ax=None, cmap='Spectral', title='$T_2^*$ (s)')[source]#

Show relaxation time distribution as stitched model section.

showWC(wlim=(0, 0.5), ax=None, cmap='Spectral', title='$\\theta$ (-)')[source]#

Show water content distribution as stitched model section.