Other functions

kernel – Kernel calculation

Kernel of empymod, calculates the wavenumber-domain electromagnetic response. Plus analytical full- and half-space solutions.

The functions wavenumber, angle_factor, fullspace, greenfct, reflections, and fields are based on source files (specified in each function) from the source code distributed with [HuTS15], which can be found at software.seg.org/2015/0001. These functions are (c) 2015 by Hunziker et al. and the Society of Exploration Geophysicists, https://software.seg.org/disclaimer.txt. Please read the NOTICE-file in the root directory for more information regarding the involved licenses.

empymod.kernel.wavenumber(zsrc, zrec, lsrc, lrec, depth, etaH, etaV, zetaH, zetaV, lambd, ab, xdirect, msrc, mrec, use_ne_eval)[source]

Calculate wavenumber domain solution.

Return the wavenumber domain solutions PJ0, PJ1, and PJ0b, which have to be transformed with a Hankel transform to the frequency domain. PJ0/PJ0b and PJ1 have to be transformed with Bessel functions of order 0 (\(J_0\)) and 1 (\(J_1\)), respectively.

This function corresponds loosely to equations 105–107, 111–116, 119–121, and 123–128 in [HuTS15], and equally loosely to the file kxwmod.c.

[HuTS15] uses Bessel functions of orders 0, 1, and 2 (\(J_0, J_1, J_2\)). The implementations of the Fast Hankel Transform and the Quadrature-with-Extrapolation in transform are set-up with Bessel functions of order 0 and 1 only. This is achieved by applying the recurrence formula

\[J_2(kr) = \frac{2}{kr} J_1(kr) - J_0(kr) \ .\]

Note

PJ0 and PJ0b could theoretically be added here into one, and then be transformed in one go. However, PJ0b has to be multiplied by factAng later. This has to be done after the Hankel transform for methods which make use of spline interpolation, in order to work for offsets that are not in line with each other.

This function is called from one of the Hankel functions in transform. Consult the modelling routines in model for a description of the input and output parameters.

If you are solely interested in the wavenumber-domain solution you can call this function directly. However, you have to make sure all input arguments are correct, as no checks are carried out here.

empymod.kernel.angle_factor(angle, ab, msrc, mrec)[source]

Return the angle-dependent factor.

The whole calculation in the wavenumber domain is only a function of the distance between the source and the receiver, it is independent of the angel. The angle-dependency is this factor, which can be applied to the corresponding parts in the wavenumber or in the frequency domain.

The angle_factor corresponds to the sine and cosine-functions in Eqs 105-107, 111-116, 119-121, 123-128.

This function is called from one of the Hankel functions in transform. Consult the modelling routines in model for a description of the input and output parameters.

empymod.kernel.fullspace(off, angle, zsrc, zrec, etaH, etaV, zetaH, zetaV, ab, msrc, mrec)[source]

Analytical full-space solutions in the frequency domain.

\[\hat{G}^{ee}_{\alpha\beta}, \hat{G}^{ee}_{3\alpha}, \hat{G}^{ee}_{33}, \hat{G}^{em}_{\alpha\beta}, \hat{G}^{em}_{\alpha 3}\]

This function corresponds to equations 45–50 in [HuTS15], and loosely to the corresponding files Gin11.F90, Gin12.F90, Gin13.F90, Gin22.F90, Gin23.F90, Gin31.F90, Gin32.F90, Gin33.F90, Gin41.F90, Gin42.F90, Gin43.F90, Gin51.F90, Gin52.F90, Gin53.F90, Gin61.F90, and Gin62.F90.

This function is called from one of the modelling routines in model. Consult these modelling routines for a description of the input and output parameters.

empymod.kernel.greenfct(zsrc, zrec, lsrc, lrec, depth, etaH, etaV, zetaH, zetaV, lambd, ab, xdirect, msrc, mrec, use_ne_eval)[source]

Calculate Green’s function for TM and TE.

\[\tilde{g}^{tm}_{hh}, \tilde{g}^{tm}_{hz}, \tilde{g}^{tm}_{zh}, \tilde{g}^{tm}_{zz}, \tilde{g}^{te}_{hh}, \tilde{g}^{te}_{zz}\]

This function corresponds to equations 108–110, 117/118, 122; 89–94, A18–A23, B13–B15; 97–102 A26–A31, and B16–B18 in [HuTS15], and loosely to the corresponding files Gamma.F90, Wprop.F90, Ptotalx.F90, Ptotalxm.F90, Ptotaly.F90, Ptotalym.F90, Ptotalz.F90, and Ptotalzm.F90.

The Green’s functions are multiplied according to Eqs 105-107, 111-116, 119-121, 123-128; with the factors inside the integrals.

This function is called from the function kernel.wavenumber.

empymod.kernel.reflections(depth, e_zH, Gam, lrec, lsrc, use_ne_eval)[source]

Calculate Rp, Rm.

\[R^\pm_n, \bar{R}^\pm_n\]

This function corresponds to equations 64/65 and A-11/A-12 in [HuTS15], and loosely to the corresponding files Rmin.F90 and Rplus.F90.

This function is called from the function kernel.greenfct.

empymod.kernel.fields(depth, Rp, Rm, Gam, lrec, lsrc, zsrc, ab, TM, use_ne_eval)[source]

Calculate Pu+, Pu-, Pd+, Pd-.

\[P^{u\pm}_s, P^{d\pm}_s, \bar{P}^{u\pm}_s, \bar{P}^{d\pm}_s; P^{u\pm}_{s-1}, P^{u\pm}_n, \bar{P}^{u\pm}_{s-1}, \bar{P}^{u\pm}_n; P^{d\pm}_{s+1}, P^{d\pm}_n, \bar{P}^{d\pm}_{s+1}, \bar{P}^{d\pm}_n\]

This function corresponds to equations 81/82, 95/96, 103/104, A-8/A-9, A-24/A-25, and A-32/A-33 in [HuTS15], and loosely to the corresponding files Pdownmin.F90, Pdownplus.F90, Pupmin.F90, and Pdownmin.F90.

This function is called from the function kernel.greenfct.

empymod.kernel.halfspace(off, angle, zsrc, zrec, etaH, etaV, freqtime, ab, signal, solution='dhs')[source]

Return frequency- or time-space domain VTI half-space solution.

Calculates the frequency- or time-space domain electromagnetic response for a half-space below air using the diffusive approximation, as given in [SlHM10], where the electric source is located at [x=0, y=0, z=zsrc>=0], and the electric receiver at [x=cos(angle)*off, y=sin(angle)*off, z=zrec>=0].

It can also be used to calculate the fullspace solution or the separate fields: direct field, reflected field, and airwave; always using the diffusive approximation. See solution-parameter.

This function is called from one of the modelling routines in model. Consult these modelling routines for a description of the input and solution parameters.

transform – Hankel and Fourier Transforms

Methods to carry out the required Hankel transform from wavenumber to frequency domain and Fourier transform from frequency to time domain.

The functions for the QWE and DLF Hankel and Fourier transforms are based on source files (specified in each function) from the source code distributed with [Key12], which can be found at software.seg.org/2012/0003. These functions are (c) 2012 by Kerry Key and the Society of Exploration Geophysicists, https://software.seg.org/disclaimer.txt. Please read the NOTICE-file in the root directory for more information regarding the involved licenses.

empymod.transform.fht(zsrc, zrec, lsrc, lrec, off, factAng, depth, ab, etaH, etaV, zetaH, zetaV, xdirect, fhtarg, use_ne_eval, msrc, mrec)[source]

Hankel Transform using the Digital Linear Filter method.

The Digital Linear Filter method was introduced to geophysics by [Ghos70], and made popular and wide-spread by [Ande75], [Ande79], [Ande82]. The DLF is sometimes referred to as the Fast Hankel Transform FHT, from which this routine has its name.

This implementation of the DLF follows [Key12], equation 6. Without going into the mathematical details (which can be found in any of the above papers) and following [Key12], the DLF method rewrites the Hankel transform of the form

\[F(r) = \int^\infty_0 f(\lambda)J_v(\lambda r)\ \mathrm{d}\lambda\]

as

\[F(r) = \sum^n_{i=1} f(b_i/r)h_i/r \ ,\]

where \(h\) is the digital filter.The Filter abscissae b is given by

\[b_i = \lambda_ir = e^{ai}, \qquad i = -l, -l+1, \cdots, l \ ,\]

with \(l=(n-1)/2\), and \(a\) is the spacing coefficient.

This function is loosely based on get_CSEM1D_FD_FHT.m from the source code distributed with [Key12].

The function is called from one of the modelling routines in model. Consult these modelling routines for a description of the input and output parameters.

Returns:
fEM : array

Returns frequency-domain EM response.

kcount : int

Kernel count. For DLF, this is 1.

conv : bool

Only relevant for QWE/QUAD.

empymod.transform.hqwe(zsrc, zrec, lsrc, lrec, off, factAng, depth, ab, etaH, etaV, zetaH, zetaV, xdirect, qweargs, use_ne_eval, msrc, mrec)[source]

Hankel Transform using Quadrature-With-Extrapolation.

Quadrature-With-Extrapolation was introduced to geophysics by [Key12]. It is one of many so-called ISE methods to solve Hankel Transforms, where ISE stands for Integration, Summation, and Extrapolation.

Following [Key12], but without going into the mathematical details here, the QWE method rewrites the Hankel transform of the form

\[F(r) = \int^\infty_0 f(\lambda)J_v(\lambda r)\ \mathrm{d}\lambda\]

as a quadrature sum which form is similar to the DLF (equation 15),

\[F_i \approx \sum^m_{j=1} f(x_j/r)w_j g(x_j) = \sum^m_{j=1} f(x_j/r)\hat{g}(x_j) \ ,\]

but with various bells and whistles applied (using the so-called Shanks transformation in the form of a routine called \(\epsilon\)-algorithm ([Shan55], [Wynn56]; implemented with algorithms from [Tref00] and [Weni89]).

This function is based on get_CSEM1D_FD_QWE.m, qwe.m, and getBesselWeights.m from the source code distributed with [Key12].

In the spline-version, hqwe checks how steep the decay of the wavenumber-domain result is, and calls QUAD for the very steep interval, for which QWE is not suited.

The function is called from one of the modelling routines in model. Consult these modelling routines for a description of the input and output parameters.

Returns:
fEM : array

Returns frequency-domain EM response.

kcount : int

Kernel count.

conv : bool

If true, QWE/QUAD converged. If not, <htarg> might have to be adjusted.

empymod.transform.hquad(zsrc, zrec, lsrc, lrec, off, factAng, depth, ab, etaH, etaV, zetaH, zetaV, xdirect, quadargs, use_ne_eval, msrc, mrec)[source]

Hankel Transform using the QUADPACK library.

This routine uses the scipy.integrate.quad module, which in turn makes use of the Fortran library QUADPACK (qagse).

It is massively (orders of magnitudes) slower than either fht or hqwe, and is mainly here for completeness and comparison purposes. It always uses interpolation in the wavenumber domain, hence it generally will not be as precise as the other methods. However, it might work in some areas where the others fail.

The function is called from one of the modelling routines in model. Consult these modelling routines for a description of the input and output parameters.

Returns:
fEM : array

Returns frequency-domain EM response.

kcount : int

Kernel count. For HQUAD, this is 1.

conv : bool

If true, QUAD converged. If not, <htarg> might have to be adjusted.

empymod.transform.ffht(fEM, time, freq, ftarg)[source]

Fourier Transform using the Digital Linear Filter method.

It follows the Filter methodology [Ande75], using Cosine- and Sine-filters; see fht for more information.

The function is called from one of the modelling routines in model. Consult these modelling routines for a description of the input and output parameters.

This function is based on get_CSEM1D_TD_FHT.m from the source code distributed with [Key12].

Returns:
tEM : array

Returns time-domain EM response of fEM for given time.

conv : bool

Only relevant for QWE/QUAD.

empymod.transform.fqwe(fEM, time, freq, qweargs)[source]

Fourier Transform using Quadrature-With-Extrapolation.

It follows the QWE methodology [Key12] for the Hankel transform, see hqwe for more information.

The function is called from one of the modelling routines in model. Consult these modelling routines for a description of the input and output parameters.

This function is based on get_CSEM1D_TD_QWE.m from the source code distributed with [Key12].

fqwe checks how steep the decay of the frequency-domain result is, and calls QUAD for the very steep interval, for which QWE is not suited.

Returns:
tEM : array

Returns time-domain EM response of fEM for given time.

conv : bool

If true, QWE/QUAD converged. If not, <ftarg> might have to be adjusted.

empymod.transform.fftlog(fEM, time, freq, ftarg)[source]

Fourier Transform using FFTLog.

FFTLog is the logarithmic analogue to the Fast Fourier Transform FFT. FFTLog was presented in Appendix B of [Hami00] and published at <http://casa.colorado.edu/~ajsh/FFTLog>.

This function uses a simplified version of pyfftlog, which is a python-version of FFTLog. For more details regarding pyfftlog see <https://github.com/prisae/pyfftlog>.

Not the full flexibility of FFTLog is available here: Only the logarithmic FFT (fftl in FFTLog), not the Hankel transform (fht in FFTLog). Furthermore, the following parameters are fixed:

  • kr = 1 (initial value)
  • kropt = 1 (silently adjusts kr)
  • dir = 1 (forward)

Furthermore, q is restricted to -1 <= q <= 1.

The function is called from one of the modelling routines in model. Consult these modelling routines for a description of the input and output parameters.

Returns:
tEM : array

Returns time-domain EM response of fEM for given time.

conv : bool

Only relevant for QWE/QUAD.

empymod.transform.fft(fEM, time, freq, ftarg)[source]

Fourier Transform using the Fast Fourier Transform.

The function is called from one of the modelling routines in model. Consult these modelling routines for a description of the input and output parameters.

Returns:
tEM : array

Returns time-domain EM response of fEM for given time.

conv : bool

Only relevant for QWE/QUAD.

empymod.transform.dlf(signal, points, out_pts, filt, pts_per_dec, kind=None, factAng=None, ab=None, int_pts=None)[source]

Digital Linear Filter method.

This is the kernel of the DLF method, used for the Hankel (fht) and the Fourier (ffht) Transforms. See fht for an extensive description.

For the Hankel transform, signal contains 3 complex wavenumber-domain signals: (PJ0, PJ1, PJ0b), as returned from kernel.wavenumber. The Hankel DLF has two additional, optional parameters: factAng, as returned from kernel.angle_factor, and ab. The PJ0-kernel is the part of the wavenumber-domain calculation which contains a zeroth-order Bessel function and does NOT depend on the angle between source and receiver, only on offset. PJ0b and PJ1 are the parts of the wavenumber-domain calculation which contain a zeroth- and first-order Bessel function, respectively, and can depend on the angle between source and receiver. PJ0, PJ1, or PJ0b can also be None, if they are not used.

For the Fourier transform, signal is a complex frequency-domain signal. The Fourier DLF requires one additional parameter, kind, which will be ‘cos’ or ‘sin’.

empymod.transform.qwe(rtol, atol, maxint, inp, intervals, lambd=None, off=None, factAng=None)[source]

Quadrature-With-Extrapolation.

This is the kernel of the QWE method, used for the Hankel (hqwe) and the Fourier (fqwe) Transforms. See hqwe for an extensive description.

This function is based on qwe.m from the source code distributed with [Key12].

empymod.transform.get_spline_values(filt, inp, nr_per_dec=None)[source]

Return required calculation points.

empymod.transform.fhti(rmin, rmax, n, q, mu)[source]

Return parameters required for FFTLog.

filters – Digital Linear Filters

Filters for the Digital Linear Filter (DLF) method for the Hankel [Ghos70]) and the Fourier ([Ande75]) transforms.

To calculate the dlf.factor I used

np.around(np.average(dlf.base[1:]/dlf.base[:-1]), 15)

The filters kong_61_2007 and kong_241_2007 from [Kong07], and key_101_2009, key_201_2009, key_401_2009, key_81_CosSin_2009, key_241_CosSin_2009, and key_601_CosSin_2009 from [Key09] are taken from DIPOLE1D, [Key09], which can be downloaded at https://marineemlab.ucsd.edu/Projects/Occam/1DCSEM (1DCSEM). DIPOLE1D is distributed under the license GNU GPL version 3 or later. Kerry Key gave his written permission to re-distribute the filters under the Apache License, Version 2.0 (email from Kerry Key to Dieter Werthmüller, 21 November 2016).

The filters anderson_801_1982 from [Ande82] and key_51_2012, key_101_2012, key_201_2012, key_101_CosSin_2012, and key_201_CosSin_2012, all from [Key12], are taken from the software distributed with [Key12] and available at https://software.seg.org/2012/0003 (SEG-2012-003). These filters are distributed under the SEG license.

The filter wer_201_2018 was designed with the add-on fdesign, see https://github.com/empymod/article-fdesign.

class empymod.filters.DigitalFilter(name, savename=None, filter_coeff=None)[source]

Simple Class for Digital Linear Filters.

Parameters:
name : str

Name of the DFL.

savename = str

Name with which the filter is saved. If None (default) it is set to the same value as name.

filter_coeff = list of str

By default, the following filter coefficients are checked:

filter_coeff = ['j0', 'j1', 'sin', 'cos']

This accounts for the standard Hankel and Fourier DLF in CSEM modelling. However, additional coefficient names can be provided via this parameter (in list format).

fromfile(self, path='filters')[source]

Load filter values from ascii-files.

Load filter base and filter coefficients from ascii files in the directory path; path can be a relative or absolute path.

Examples

>>> import empymod
>>> # Create an empty filter;
>>> # Name has to be the base of the text files
>>> filt = empymod.filters.DigitalFilter('my-filter')
>>> # Load the ascii-files
>>> filt.fromfile()
>>> # This will load the following three files:
>>> #    ./filters/my-filter_base.txt
>>> #    ./filters/my-filter_j0.txt
>>> #    ./filters/my-filter_j1.txt
>>> # and store them in filt.base, filt.j0, and filt.j1.
tofile(self, path='filters')[source]

Save filter values to ascii-files.

Store the filter base and the filter coefficients in separate files in the directory path; path can be a relative or absolute path.

Examples

>>> import empymod
>>> # Load a filter
>>> filt = empymod.filters.wer_201_2018()
>>> # Save it to pure ascii-files
>>> filt.tofile()
>>> # This will save the following three files:
>>> #    ./filters/wer_201_2018_base.txt
>>> #    ./filters/wer_201_2018_j0.txt
>>> #    ./filters/wer_201_2018_j1.txt
empymod.filters.kong_61_2007()[source]

Kong 61 pt Hankel filter, as published in [Kong07].

Taken from file FilterModules.f90 provided with 1DCSEM.

License: Apache License, Version 2.0,.

empymod.filters.kong_241_2007()[source]

Kong 241 pt Hankel filter, as published in [Kong07].

Taken from file FilterModules.f90 provided with 1DCSEM.

License: Apache License, Version 2.0,.

empymod.filters.key_101_2009()[source]

Key 101 pt Hankel filter, as published in [Key09].

Taken from file FilterModules.f90 provided with 1DCSEM.

License: Apache License, Version 2.0,.

empymod.filters.key_201_2009()[source]

Key 201 pt Hankel filter, as published in [Key09].

Taken from file FilterModules.f90 provided with 1DCSEM.

License: Apache License, Version 2.0,.

empymod.filters.key_401_2009()[source]

Key 401 pt Hankel filter, as published in [Key09].

Taken from file FilterModules.f90 provided with 1DCSEM.

License: Apache License, Version 2.0,.

empymod.filters.anderson_801_1982()[source]

Anderson 801 pt Hankel filter, as published in [Ande82].

Taken from file wa801Hankel.txt provided with SEG-2012-003.

License: https://software.seg.org/disclaimer.txt.

empymod.filters.key_51_2012()[source]

Key 51 pt Hankel filter, as published in [Key12].

Taken from file kk51Hankel.txt provided with SEG-2012-003.

License: https://software.seg.org/disclaimer.txt.

empymod.filters.key_101_2012()[source]

Key 101 pt Hankel filter, as published in [Key12].

Taken from file kk101Hankel.txt provided with SEG-2012-003.

License: https://software.seg.org/disclaimer.txt.

empymod.filters.key_201_2012()[source]

Key 201 pt Hankel filter, as published in [Key12].

Taken from file kk201Hankel.txt provided with SEG-2012-003.

License: https://software.seg.org/disclaimer.txt.

empymod.filters.wer_201_2018()[source]

Werthmüller 201 pt Hankel filter, 2018.

Designed with the empymod add-on fdesign, see https://github.com/empymod/article-fdesign.

License: Apache License, Version 2.0,.

empymod.filters.key_81_CosSin_2009()[source]

Key 81 pt CosSin filter, as published in [Key09].

Taken from file FilterModules.f90 provided with 1DCSEM.

License: Apache License, Version 2.0,.

empymod.filters.key_241_CosSin_2009()[source]

Key 241 pt CosSin filter, as published in [Key09].

Taken from file FilterModules.f90 provided with 1DCSEM.

License: Apache License, Version 2.0,.

empymod.filters.key_601_CosSin_2009()[source]

Key 601 pt CosSin filter, as published in [Key09].

Taken from file FilterModules.f90 provided with 1DCSEM.

License: Apache License, Version 2.0,.

empymod.filters.key_101_CosSin_2012()[source]

Key 101 pt CosSin filter, as published in [Key12].

Taken from file kk101CosSin.txt provided with SEG-2012-003.

License: https://software.seg.org/disclaimer.txt.

empymod.filters.key_201_CosSin_2012()[source]

Key 201 pt CosSin filter, as published in [Key12].

Taken from file kk201CosSin.txt provided with SEG-2012-003.

License: https://software.seg.org/disclaimer.txt.

utils – Utilites

Utilities for model such as checking input parameters.

This module consists of four groups of functions:
  1. General settings
  2. Class EMArray
  3. Input parameter checks for modelling
  4. Internal utilities
class empymod.utils.EMArray[source]

Subclassing an ndarray: add amplitude <amp> and phase <pha>.

Parameters:
data : array

Data to which to add .amp and .pha attributes.

Examples

>>> import numpy as np
>>> from empymod.utils import EMArray
>>> emvalues = EMArray(np.array([1,2,3])+1j*np.array([1, 0, -1]))
>>> print('Amplitude : ', emvalues.amp)
Amplitude :  [ 1.41421356  2.          3.16227766]
>>> print('Phase     : ', emvalues.pha)
Phase     :  [ 45.           0.         -18.43494882]
Attributes:
amp : ndarray

Make amplitude an attribute.

pha : ndarray

Make phase an attribute (unwrapped and in degrees.

amp

Make amplitude an attribute.

pha

Make phase an attribute (unwrapped and in degrees.

empymod.utils.check_time_only(time, signal, verb)[source]

Check time and signal parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
time : array_like

Times t (s).

signal : {None, 0, 1, -1}
Source signal:
  • None: Frequency-domain response
  • -1 : Switch-off time-domain response
  • 0 : Impulse time-domain response
  • +1 : Switch-on time-domain response
verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
time : float

Time, checked for size and assured min_time.

empymod.utils.check_time(time, signal, ft, ftarg, verb)[source]

Check time domain specific input parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
time : array_like

Times t (s).

signal : {None, 0, 1, -1}
Source signal:
  • None: Frequency-domain response
  • -1 : Switch-off time-domain response
  • 0 : Impulse time-domain response
  • +1 : Switch-on time-domain response
ft : {‘sin’, ‘cos’, ‘qwe’, ‘fftlog’, ‘fft’}

Flag for Fourier transform.

ftarg : str or filter from empymod.filters or array_like,

Only used if signal !=None. Depends on the value for ft:

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
time : float

Time, checked for size and assured min_time.

freq : float

Frequencies required for given times and ft-settings.

ft, ftarg

Checked if valid and set to defaults if not provided, checked with signal.

empymod.utils.check_model(depth, res, aniso, epermH, epermV, mpermH, mpermV, xdirect, verb)[source]

Check the model: depth and corresponding layer parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
depth : list

Absolute layer interfaces z (m); #depth = #res - 1 (excluding +/- infinity).

res : array_like

Horizontal resistivities rho_h (Ohm.m); #res = #depth + 1.

aniso : array_like

Anisotropies lambda = sqrt(rho_v/rho_h) (-); #aniso = #res.

epermH, epermV : array_like

Relative horizontal/vertical electric permittivities epsilon_h/epsilon_v (-); #epermH = #epermV = #res.

mpermH, mpermV : array_like

Relative horizontal/vertical magnetic permeabilities mu_h/mu_v (-); #mpermH = #mpermV = #res.

xdirect : bool, optional

If True and source and receiver are in the same layer, the direct field is calculated analytically in the frequency domain, if False it is calculated in the wavenumber domain.

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
depth : array

Depths of layer interfaces, adds -infty at beginning if not present.

res : array

As input, checked for size.

aniso : array

As input, checked for size. If None, defaults to an array of ones.

epermH, epermV : array_like

As input, checked for size. If None, defaults to an array of ones.

mpermH, mpermV : array_like

As input, checked for size. If None, defaults to an array of ones.

isfullspace : bool

If True, the model is a fullspace (res, aniso, epermH, epermV, mpermM, and mpermV are in all layers the same).

empymod.utils.check_frequency(freq, res, aniso, epermH, epermV, mpermH, mpermV, verb)[source]

Calculate frequency-dependent parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
freq : array_like

Frequencies f (Hz).

res : array_like

Horizontal resistivities rho_h (Ohm.m); #res = #depth + 1.

aniso : array_like

Anisotropies lambda = sqrt(rho_v/rho_h) (-); #aniso = #res.

epermH, epermV : array_like

Relative horizontal/vertical electric permittivities epsilon_h/epsilon_v (-); #epermH = #epermV = #res.

mpermH, mpermV : array_like

Relative horizontal/vertical magnetic permeabilities mu_h/mu_v (-); #mpermH = #mpermV = #res.

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
freq : float

Frequency, checked for size and assured min_freq.

etaH, etaV : array

Parameters etaH/etaV, same size as provided resistivity.

zetaH, zetaV : array

Parameters zetaH/zetaV, same size as provided resistivity.

empymod.utils.check_hankel(ht, htarg, verb)[source]

Check Hankel transform parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
ht : {‘fht’, ‘qwe’, ‘quad’}

Flag to choose the Hankel transform.

htarg : str or filter from empymod.filters or array_like,

Depends on the value for ht.

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
ht, htarg

Checked if valid and set to defaults if not provided.

empymod.utils.check_opt(opt, loop, ht, htarg, verb)[source]

Check optimization parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
opt : {None, ‘parallel’}

Optimization flag; use numexpr or not.

loop : {None, ‘freq’, ‘off’}

Loop flag.

ht : str

Flag to choose the Hankel transform.

htarg : array_like,

Depends on the value for ht.

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
use_ne_eval : bool

Boolean if to use numexpr.

loop_freq : bool

Boolean if to loop over frequencies.

loop_off : bool

Boolean if to loop over offsets.

empymod.utils.check_dipole(inp, name, verb)[source]

Check dipole parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
inp : list of floats or arrays

Pole coordinates (m): [pole-x, pole-y, pole-z].

name : str, {‘src’, ‘rec’}

Pole-type.

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
inp : list

List of pole coordinates [x, y, z].

ninp : int

Number of inp-elements

empymod.utils.check_bipole(inp, name)[source]

Check di-/bipole parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
inp : list of floats or arrays

Coordinates of inp (m): [dipole-x, dipole-y, dipole-z, azimuth, dip] or. [bipole-x0, bipole-x1, bipole-y0, bipole-y1, bipole-z0, bipole-z1].

name : str, {‘src’, ‘rec’}

Pole-type.

Returns:
inp : list

As input, checked for type and length.

ninp : int

Number of inp.

ninpz : int

Number of inp depths (ninpz is either 1 or ninp).

isdipole : bool

True if inp is a dipole.

empymod.utils.check_ab(ab, verb)[source]

Check source-receiver configuration.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
ab : int

Source-receiver configuration.

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
ab_calc : int

Adjusted source-receiver configuration using reciprocity.

msrc, mrec : bool

If True, src/rec is magnetic; if False, src/rec is electric.

empymod.utils.check_solution(solution, signal, ab, msrc, mrec)[source]

Check required solution with parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
solution : str

String to define analytical solution.

signal : {None, 0, 1, -1}
Source signal:
  • None: Frequency-domain response
  • -1 : Switch-off time-domain response
  • 0 : Impulse time-domain response
  • +1 : Switch-on time-domain response
msrc, mrec : bool

True if src/rec is magnetic, else False.

empymod.utils.get_abs(msrc, mrec, srcazm, srcdip, recazm, recdip, verb)[source]

Get required ab’s for given angles.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
msrc, mrec : bool

True if src/rec is magnetic, else False.

srcazm, recazm : float

Horizontal source/receiver angle (azimuth).

srcdip, recdip : float

Vertical source/receiver angle (dip).

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
ab_calc : array of int

ab’s to calculate for this bipole.

empymod.utils.get_geo_fact(ab, srcazm, srcdip, recazm, recdip, msrc, mrec)[source]

Get required geometrical scaling factor for given angles.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
ab : int

Source-receiver configuration.

srcazm, recazm : float

Horizontal source/receiver angle.

srcdip, recdip : float

Vertical source/receiver angle.

Returns:
fact : float

Geometrical scaling factor.

empymod.utils.get_azm_dip(inp, iz, ninpz, intpts, isdipole, strength, name, verb)[source]

Get angles, interpolation weights and normalization weights.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
inp : list of floats or arrays
Input coordinates (m):
  • [x0, x1, y0, y1, z0, z1] (bipole of finite length)
  • [x, y, z, azimuth, dip] (dipole, infinitesimal small)
iz : int

Index of current di-/bipole depth (-).

ninpz : int

Total number of di-/bipole depths (ninpz = 1 or npinz = nsrc) (-).

intpts : int

Number of integration points for bipole (-).

isdipole : bool

Boolean if inp is a dipole.

strength : float, optional
Source strength (A):
  • If 0, output is normalized to source and receiver of 1 m length, and source strength of 1 A.
  • If != 0, output is returned for given source and receiver length, and source strength.
name : str, {‘src’, ‘rec’}

Pole-type.

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
tout : list of floats or arrays

Dipole coordinates x, y, and z (m).

azm : float or array of floats

Horizontal angle (azimuth).

dip : float or array of floats

Vertical angle (dip).

g_w : float or array of floats

Factors from Gaussian interpolation.

intpts : int

As input, checked.

inp_w : float or array of floats

Factors from source/receiver length and source strength.

empymod.utils.get_off_ang(src, rec, nsrc, nrec, verb)[source]

Get depths, offsets, angles, hence spatial input parameters.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
src, rec : list of floats or arrays

Source/receiver dipole coordinates x, y, and z (m).

nsrc, nrec : int

Number of sources/receivers (-).

verb : {0, 1, 2, 3, 4}

Level of verbosity.

Returns:
off : array of floats

Offsets

angle : array of floats

Angles

empymod.utils.get_layer_nr(inp, depth)[source]

Get number of layer in which inp resides.

Note: If zinp is on a layer interface, the layer above the interface is chosen.

This check-function is called from one of the modelling routines in model. Consult these modelling routines for a detailed description of the input parameters.

Parameters:
inp : list of floats or arrays

Dipole coordinates (m)

depth : array

Depths of layer interfaces.

Returns:
linp : int or array_like of int

Layer number(s) in which inp resides (plural only if bipole).

zinp : float or array

inp[2] (depths).

empymod.utils.printstartfinish(verb, inp=None, kcount=None)[source]

Print start and finish with time measure and kernel count.

empymod.utils.conv_warning(conv, targ, name, verb)[source]

Print error if QWE/QUAD did not converge at least once.

empymod.utils.set_minimum(min_freq=None, min_time=None, min_off=None, min_res=None, min_angle=None)[source]

Set minimum values of parameters.

The given parameters are set to its minimum value if they are smaller.

Parameters:
min_freq : float, optional

Minimum frequency [Hz] (default 1e-20 Hz).

min_time : float, optional

Minimum time [s] (default 1e-20 s).

min_off : float, optional

Minimum offset [m] (default 1e-3 m). Also used to round src- & rec-coordinates.

min_res : float, optional

Minimum horizontal and vertical resistivity [Ohm.m] (default 1e-20).

min_angle : float, optional

Minimum angle [-] (default 1e-10).

empymod.utils.get_minimum()[source]

Return the current minimum values.

Returns:
min_vals : dict

Dictionary of current minimum values with keys

  • min_freq : float
  • min_time : float
  • min_off : float
  • min_res : float
  • min_angle : float

For a full description of these options, see set_minimum.

empymod.utils.spline_backwards_hankel(ht, htarg, opt)[source]

Check opt if deprecated ‘spline’ is used.

Returns corrected htarg, opt. r

class empymod.utils.Report(add_pckg=None, ncol=3, text_width=80, sort=False)[source]

Print date, time, and version information.

Use scooby to print date, time, and package version information in any environment (Jupyter notebook, IPython console, Python console, QT console), either as html-table (notebook) or as plain text (anywhere).

Always shown are the OS, number of CPU(s), numpy, scipy, empymod, sys.version, and time/date.

Additionally shown are, if they can be imported, numexpr, IPython, and matplotlib. It also shows MKL information, if available.

All modules provided in add_pckg are also shown.

Parameters:
add_pckg : packages, optional

Package or list of packages to add to output information (must be imported beforehand).

ncol : int, optional

Number of package-columns in html table (no effect in text-version); Defaults to 3.

text_width : int, optional

The text width for non-HTML display modes

sort : bool, optional

Sort the packages when the report is shown

Examples

>>> import pytest
>>> import dateutil
>>> from emg3d import Report
>>> Report()                            # Default values
>>> Report(pytest)                      # Provide additional package
>>> Report([pytest, dateutil], ncol=5)  # Set nr of columns
class empymod.utils.Versions(add_pckg=None, ncol=3)[source]

New name is Report, here for backwards compatibility.

empymod.utils.versions(mode=None, add_pckg=None, ncol=4)[source]

Old func-way of class Report, here for backwards compatibility.