Equivalent layer processing (fatiando.gravmag.eqlayer)

Equivalent layer processing.

Use the classes here to estimate an equivalent layer from potential field data. Then you can use the estimated layer to perform tranformations (gridding, continuation, derivation, reduction to the pole, etc.) by forward modeling the layer. Use fatiando.gravmag.sphere for forward modeling.

Algorithms

  • EQLGravity and EQLTotalField: The classic (space domain) equivalent layer as formulated in Li and Oldenburg (2010) or Oliveira Jr. et al (2012). Doesn’t have wavelet compression or other tweaks.
  • PELGravity and PELTotalField: The polynomial equivalent layer of Oliveira Jr. et al (2012). A fast and memory efficient algorithm. Both of these require special regularization (PELSmoothness).

References

Li, Y., and D. W. Oldenburg (2010), Rapid construction of equivalent sources using wavelets, Geophysics, 75(3), L51-L59, doi:10.1190/1.3378764.

Oliveira Jr., V. C., V. C. F. Barbosa, and L. Uieda (2012), Polynomial equivalent layer, Geophysics, 78(1), G1-G13, doi:10.1190/geo2012-0196.1.


class fatiando.gravmag.eqlayer.EQLBase(x, y, z, data, grid)[source]

Bases: fatiando.inversion.misfit.Misfit

Base class for the classic equivalent layer.

config(method, **kwargs)

Configure the optimization method and its parameters.

This sets the method used by fit and the keyword arguments that are passed to it.

Parameters:

  • method
    : string

    The optimization method. One of: 'linear', 'newton', 'levmarq', 'steepest', 'acor'

Other keyword arguments that can be passed are the ones allowed by each method.

Some methods have required arguments:

  • newton, levmarq and steepest require the initial argument (an initial estimate for the gradient descent)
  • acor requires the bounds argument (min/max values for the search space)

See the corresponding docstrings for more information:

copy(deep=False)

Make a copy of me together with all the cached methods.

estimate_

A nicely formatted version of the estimate.

If the class implements a fmt_estimate method, this will its results. This can be used to convert the parameter vector to a more useful form, like a fatiando.mesher object.

fit()

Solve for the parameter vector that minimizes this objective function.

Uses the optimization method and parameters defined using the config method.

The estimated parameter vector can be accessed through the p_ attribute. A (possibly) formatted version (converted to a more manageable type) of the estimate can be accessed through the property estimate_.

fmt_estimate(p)

Called when accessing the property estimate_.

Use this to convert the parameter vector (p) to a more useful form, like a geometric object, etc.

Parameters:

  • p
    : 1d-array

    The parameter vector.

Returns:

  • formatted

    Pretty much anything you want.

gradient(p)

The gradient vector of the misfit function.

\[\bar{g} = -2\bar{\bar{J}}^T\bar{r}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix and \(\bar{r}\) is the residual vector.

Parameters:

  • p
    : 1d-array

    The parameter vector where the gradient is evaluated

Returns:

  • gradient
    : 1d-array

    The gradient vector.

hessian(p)

The Hessian of the misfit function with respect to the parameters.

Calculated using the Gauss approximation:

\[\bar{\bar{H}} \approx 2\bar{\bar{J}}^T\bar{\bar{J}}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix.

For linear problems, the Hessian matrix is cached in memory, so calling this method again will not trigger a re-calculation.

Parameters:

  • p
    : 1d-array

    The parameter vector where the Hessian is evaluated

Returns:

  • hessian
    : 2d-array

    The Hessian matrix

predicted(p)[source]

Calculate the data predicted by a given parameter vector.

Parameters:

  • p
    : 1d-array (optional)

    The parameter vector with the estimated physical properties of the layer. If not given, will use the value calculated by .fit().

Returns:

  • result
    : 1d-array

    The predicted data vector.

regul_param

The regularization parameter (scale factor) for the objetive function.

Defaults to 1.

residuals(p=None)

Calculate the residuals vector (observed - predicted data).

Parameters:

  • p
    : 1d-array or None

    The parameter vector used to calculate the residuals. If None, will use the current estimate stored in estimate_.

Returns:

  • residuals
    : 1d-array or list of 1d-arrays

    The residual vector. If this is the sum of 1 or more Misfit instances, will return the residual vector from each of the summed misfits in the order of the sum.

set_weights(weights)

Set the data weights.

Using weights for the data, the least-squares data-misfit function becomes:

\[\phi = \bar{r}^T \bar{\bar{W}}\bar{r}\]

Parameters:

  • weights
    : 1d-array or 2d-array or None

    Weights for the data vector. If None, will remove any weights that have been set before. If it is a 2d-array, it will be interpreted as the weight matrix \(\bar{\bar{W}}\). If it is a 1d-array, it will be interpreted as the diagonal of the weight matrix (all off-diagonal elements will default to zero). The weight matrix can be a sparse array from scipy.sparse.

value(p)

Calculate the value of the misfit for a given parameter vector.

The value is given by:

\[\phi(\bar{p}) = \bar{r}^T\bar{\bar{W}}\bar{r}\]

where \(\bar{r}\) is the residual vector and \(bar{\bar{W}}\) are optional data weights.

Parameters:

  • p
    : 1d-array or None

    The parameter vector.

Returns:

  • value
    : float

    The value of the misfit function.

class fatiando.gravmag.eqlayer.EQLGravity(x, y, z, data, grid, field='gz')[source]

Bases: fatiando.gravmag.eqlayer.EQLBase

Estimate an equivalent layer from gravity data.

Note

Assumes x = North, y = East, z = Down.

Parameters:

  • x, y, z
    : 1d-arrays

    The x, y, z coordinates of each data point.

  • data
    : 1d-array

    The gravity data at each point.

  • grid
    : PointGrid

    The sources in the equivalent layer. Will invert for the density of each point in the grid.

  • field
    : string

    Which gravitational field is the data. Options are: 'gz' (gravity anomaly), 'gxx', 'gxy', ..., 'gzz' (gravity gradient tensor). Defaults to 'gz'.

config(method, **kwargs)

Configure the optimization method and its parameters.

This sets the method used by fit and the keyword arguments that are passed to it.

Parameters:

  • method
    : string

    The optimization method. One of: 'linear', 'newton', 'levmarq', 'steepest', 'acor'

Other keyword arguments that can be passed are the ones allowed by each method.

Some methods have required arguments:

  • newton, levmarq and steepest require the initial argument (an initial estimate for the gradient descent)
  • acor requires the bounds argument (min/max values for the search space)

See the corresponding docstrings for more information:

copy(deep=False)

Make a copy of me together with all the cached methods.

estimate_

A nicely formatted version of the estimate.

If the class implements a fmt_estimate method, this will its results. This can be used to convert the parameter vector to a more useful form, like a fatiando.mesher object.

fit()

Solve for the parameter vector that minimizes this objective function.

Uses the optimization method and parameters defined using the config method.

The estimated parameter vector can be accessed through the p_ attribute. A (possibly) formatted version (converted to a more manageable type) of the estimate can be accessed through the property estimate_.

fmt_estimate(p)

Called when accessing the property estimate_.

Use this to convert the parameter vector (p) to a more useful form, like a geometric object, etc.

Parameters:

  • p
    : 1d-array

    The parameter vector.

Returns:

  • formatted

    Pretty much anything you want.

gradient(p)

The gradient vector of the misfit function.

\[\bar{g} = -2\bar{\bar{J}}^T\bar{r}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix and \(\bar{r}\) is the residual vector.

Parameters:

  • p
    : 1d-array

    The parameter vector where the gradient is evaluated

Returns:

  • gradient
    : 1d-array

    The gradient vector.

hessian(p)

The Hessian of the misfit function with respect to the parameters.

Calculated using the Gauss approximation:

\[\bar{\bar{H}} \approx 2\bar{\bar{J}}^T\bar{\bar{J}}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix.

For linear problems, the Hessian matrix is cached in memory, so calling this method again will not trigger a re-calculation.

Parameters:

  • p
    : 1d-array

    The parameter vector where the Hessian is evaluated

Returns:

  • hessian
    : 2d-array

    The Hessian matrix

jacobian(p)[source]

Calculate the Jacobian matrix for a given parameter vector.

predicted(p)

Calculate the data predicted by a given parameter vector.

Parameters:

  • p
    : 1d-array (optional)

    The parameter vector with the estimated physical properties of the layer. If not given, will use the value calculated by .fit().

Returns:

  • result
    : 1d-array

    The predicted data vector.

regul_param

The regularization parameter (scale factor) for the objetive function.

Defaults to 1.

residuals(p=None)

Calculate the residuals vector (observed - predicted data).

Parameters:

  • p
    : 1d-array or None

    The parameter vector used to calculate the residuals. If None, will use the current estimate stored in estimate_.

Returns:

  • residuals
    : 1d-array or list of 1d-arrays

    The residual vector. If this is the sum of 1 or more Misfit instances, will return the residual vector from each of the summed misfits in the order of the sum.

set_weights(weights)

Set the data weights.

Using weights for the data, the least-squares data-misfit function becomes:

\[\phi = \bar{r}^T \bar{\bar{W}}\bar{r}\]

Parameters:

  • weights
    : 1d-array or 2d-array or None

    Weights for the data vector. If None, will remove any weights that have been set before. If it is a 2d-array, it will be interpreted as the weight matrix \(\bar{\bar{W}}\). If it is a 1d-array, it will be interpreted as the diagonal of the weight matrix (all off-diagonal elements will default to zero). The weight matrix can be a sparse array from scipy.sparse.

value(p)

Calculate the value of the misfit for a given parameter vector.

The value is given by:

\[\phi(\bar{p}) = \bar{r}^T\bar{\bar{W}}\bar{r}\]

where \(\bar{r}\) is the residual vector and \(bar{\bar{W}}\) are optional data weights.

Parameters:

  • p
    : 1d-array or None

    The parameter vector.

Returns:

  • value
    : float

    The value of the misfit function.

class fatiando.gravmag.eqlayer.EQLTotalField(x, y, z, data, inc, dec, grid, sinc=None, sdec=None)[source]

Bases: fatiando.gravmag.eqlayer.EQLBase

Estimate an equivalent layer from total field magnetic anomaly data.

Note

Assumes x = North, y = East, z = Down.

Parameters:

  • x, y, z
    : 1d-arrays

    The x, y, z coordinates of each data point.

  • data
    : 1d-array

    The total field anomaly data at each point.

  • inc, dec
    : floats

    The inclination and declination of the inducing field

  • grid
    : PointGrid

    The sources in the equivalent layer. Will invert for the magnetization intensity of each point in the grid.

  • sinc, sdec
    : None or floats

    The inclination and declination of the equivalent layer. Use these if there is remanent magnetization and the total magnetization of the layer if different from the induced magnetization. If there is only induced magnetization, use None

config(method, **kwargs)

Configure the optimization method and its parameters.

This sets the method used by fit and the keyword arguments that are passed to it.

Parameters:

  • method
    : string

    The optimization method. One of: 'linear', 'newton', 'levmarq', 'steepest', 'acor'

Other keyword arguments that can be passed are the ones allowed by each method.

Some methods have required arguments:

  • newton, levmarq and steepest require the initial argument (an initial estimate for the gradient descent)
  • acor requires the bounds argument (min/max values for the search space)

See the corresponding docstrings for more information:

copy(deep=False)

Make a copy of me together with all the cached methods.

estimate_

A nicely formatted version of the estimate.

If the class implements a fmt_estimate method, this will its results. This can be used to convert the parameter vector to a more useful form, like a fatiando.mesher object.

fit()

Solve for the parameter vector that minimizes this objective function.

Uses the optimization method and parameters defined using the config method.

The estimated parameter vector can be accessed through the p_ attribute. A (possibly) formatted version (converted to a more manageable type) of the estimate can be accessed through the property estimate_.

fmt_estimate(p)

Called when accessing the property estimate_.

Use this to convert the parameter vector (p) to a more useful form, like a geometric object, etc.

Parameters:

  • p
    : 1d-array

    The parameter vector.

Returns:

  • formatted

    Pretty much anything you want.

gradient(p)

The gradient vector of the misfit function.

\[\bar{g} = -2\bar{\bar{J}}^T\bar{r}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix and \(\bar{r}\) is the residual vector.

Parameters:

  • p
    : 1d-array

    The parameter vector where the gradient is evaluated

Returns:

  • gradient
    : 1d-array

    The gradient vector.

hessian(p)

The Hessian of the misfit function with respect to the parameters.

Calculated using the Gauss approximation:

\[\bar{\bar{H}} \approx 2\bar{\bar{J}}^T\bar{\bar{J}}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix.

For linear problems, the Hessian matrix is cached in memory, so calling this method again will not trigger a re-calculation.

Parameters:

  • p
    : 1d-array

    The parameter vector where the Hessian is evaluated

Returns:

  • hessian
    : 2d-array

    The Hessian matrix

jacobian(p)[source]

Calculate the Jacobian matrix for a given parameter vector.

predicted(p)

Calculate the data predicted by a given parameter vector.

Parameters:

  • p
    : 1d-array (optional)

    The parameter vector with the estimated physical properties of the layer. If not given, will use the value calculated by .fit().

Returns:

  • result
    : 1d-array

    The predicted data vector.

regul_param

The regularization parameter (scale factor) for the objetive function.

Defaults to 1.

residuals(p=None)

Calculate the residuals vector (observed - predicted data).

Parameters:

  • p
    : 1d-array or None

    The parameter vector used to calculate the residuals. If None, will use the current estimate stored in estimate_.

Returns:

  • residuals
    : 1d-array or list of 1d-arrays

    The residual vector. If this is the sum of 1 or more Misfit instances, will return the residual vector from each of the summed misfits in the order of the sum.

set_weights(weights)

Set the data weights.

Using weights for the data, the least-squares data-misfit function becomes:

\[\phi = \bar{r}^T \bar{\bar{W}}\bar{r}\]

Parameters:

  • weights
    : 1d-array or 2d-array or None

    Weights for the data vector. If None, will remove any weights that have been set before. If it is a 2d-array, it will be interpreted as the weight matrix \(\bar{\bar{W}}\). If it is a 1d-array, it will be interpreted as the diagonal of the weight matrix (all off-diagonal elements will default to zero). The weight matrix can be a sparse array from scipy.sparse.

value(p)

Calculate the value of the misfit for a given parameter vector.

The value is given by:

\[\phi(\bar{p}) = \bar{r}^T\bar{\bar{W}}\bar{r}\]

where \(\bar{r}\) is the residual vector and \(bar{\bar{W}}\) are optional data weights.

Parameters:

  • p
    : 1d-array or None

    The parameter vector.

Returns:

  • value
    : float

    The value of the misfit function.

class fatiando.gravmag.eqlayer.PELBase(x, y, z, data, grid, windows, degree)[source]

Bases: fatiando.gravmag.eqlayer.EQLBase

Base class for the Polynomial Equivalent Layer.

Note

Overloads fit to convert the estimated coefficients to physical properties. The coefficients are stored in the coeffs_ attribute.

config(method, **kwargs)

Configure the optimization method and its parameters.

This sets the method used by fit and the keyword arguments that are passed to it.

Parameters:

  • method
    : string

    The optimization method. One of: 'linear', 'newton', 'levmarq', 'steepest', 'acor'

Other keyword arguments that can be passed are the ones allowed by each method.

Some methods have required arguments:

  • newton, levmarq and steepest require the initial argument (an initial estimate for the gradient descent)
  • acor requires the bounds argument (min/max values for the search space)

See the corresponding docstrings for more information:

copy(deep=False)

Make a copy of me together with all the cached methods.

estimate_

A nicely formatted version of the estimate.

If the class implements a fmt_estimate method, this will its results. This can be used to convert the parameter vector to a more useful form, like a fatiando.mesher object.

fit()

Solve for the parameter vector that minimizes this objective function.

Uses the optimization method and parameters defined using the config method.

The estimated parameter vector can be accessed through the p_ attribute. A (possibly) formatted version (converted to a more manageable type) of the estimate can be accessed through the property estimate_.

fmt_estimate(coefs)[source]

Convert the estimated polynomial coefficients to physical property values along the layer.

Parameters:

  • coefs
    : 1d-array

    The estimated parameter vector with the polynomial coefficients

Returns:

  • estimate
    : 1d-array

    The converted physical property values along the layer.

gradient(p)

The gradient vector of the misfit function.

\[\bar{g} = -2\bar{\bar{J}}^T\bar{r}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix and \(\bar{r}\) is the residual vector.

Parameters:

  • p
    : 1d-array

    The parameter vector where the gradient is evaluated

Returns:

  • gradient
    : 1d-array

    The gradient vector.

hessian(p)

The Hessian of the misfit function with respect to the parameters.

Calculated using the Gauss approximation:

\[\bar{\bar{H}} \approx 2\bar{\bar{J}}^T\bar{\bar{J}}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix.

For linear problems, the Hessian matrix is cached in memory, so calling this method again will not trigger a re-calculation.

Parameters:

  • p
    : 1d-array

    The parameter vector where the Hessian is evaluated

Returns:

  • hessian
    : 2d-array

    The Hessian matrix

predicted(p)

Calculate the data predicted by a given parameter vector.

Parameters:

  • p
    : 1d-array (optional)

    The parameter vector with the estimated physical properties of the layer. If not given, will use the value calculated by .fit().

Returns:

  • result
    : 1d-array

    The predicted data vector.

regul_param

The regularization parameter (scale factor) for the objetive function.

Defaults to 1.

residuals(p=None)

Calculate the residuals vector (observed - predicted data).

Parameters:

  • p
    : 1d-array or None

    The parameter vector used to calculate the residuals. If None, will use the current estimate stored in estimate_.

Returns:

  • residuals
    : 1d-array or list of 1d-arrays

    The residual vector. If this is the sum of 1 or more Misfit instances, will return the residual vector from each of the summed misfits in the order of the sum.

set_weights(weights)

Set the data weights.

Using weights for the data, the least-squares data-misfit function becomes:

\[\phi = \bar{r}^T \bar{\bar{W}}\bar{r}\]

Parameters:

  • weights
    : 1d-array or 2d-array or None

    Weights for the data vector. If None, will remove any weights that have been set before. If it is a 2d-array, it will be interpreted as the weight matrix \(\bar{\bar{W}}\). If it is a 1d-array, it will be interpreted as the diagonal of the weight matrix (all off-diagonal elements will default to zero). The weight matrix can be a sparse array from scipy.sparse.

value(p)

Calculate the value of the misfit for a given parameter vector.

The value is given by:

\[\phi(\bar{p}) = \bar{r}^T\bar{\bar{W}}\bar{r}\]

where \(\bar{r}\) is the residual vector and \(bar{\bar{W}}\) are optional data weights.

Parameters:

  • p
    : 1d-array or None

    The parameter vector.

Returns:

  • value
    : float

    The value of the misfit function.

class fatiando.gravmag.eqlayer.PELGravity(x, y, z, data, grid, windows, degree, field='gz')[source]

Bases: fatiando.gravmag.eqlayer.PELBase

Estimate a polynomial equivalent layer from gravity data.

Note

Assumes x = North, y = East, z = Down.

Parameters:

  • x, y, z
    : 1d-arrays

    The x, y, z coordinates of each data point.

  • data
    : 1d-array

    The gravity data at each point.

  • grid
    : PointGrid

    The sources in the equivalent layer. Will invert for the density of each point in the grid.

  • windows
    : tuple = (ny, nx)

    The number of windows that the layer will be divided in the y and x directions, respectively

  • degree
    : int

    The degree of the bivariate polynomials used in each window of the PEL

  • field
    : string

    Which gravitational field is the data. Options are: 'gz' (gravity anomaly), 'gxx', 'gxy', ..., 'gzz' (gravity gradient tensor). Defaults to 'gz'.

config(method, **kwargs)

Configure the optimization method and its parameters.

This sets the method used by fit and the keyword arguments that are passed to it.

Parameters:

  • method
    : string

    The optimization method. One of: 'linear', 'newton', 'levmarq', 'steepest', 'acor'

Other keyword arguments that can be passed are the ones allowed by each method.

Some methods have required arguments:

  • newton, levmarq and steepest require the initial argument (an initial estimate for the gradient descent)
  • acor requires the bounds argument (min/max values for the search space)

See the corresponding docstrings for more information:

copy(deep=False)

Make a copy of me together with all the cached methods.

estimate_

A nicely formatted version of the estimate.

If the class implements a fmt_estimate method, this will its results. This can be used to convert the parameter vector to a more useful form, like a fatiando.mesher object.

fit()

Solve for the parameter vector that minimizes this objective function.

Uses the optimization method and parameters defined using the config method.

The estimated parameter vector can be accessed through the p_ attribute. A (possibly) formatted version (converted to a more manageable type) of the estimate can be accessed through the property estimate_.

fmt_estimate(coefs)

Convert the estimated polynomial coefficients to physical property values along the layer.

Parameters:

  • coefs
    : 1d-array

    The estimated parameter vector with the polynomial coefficients

Returns:

  • estimate
    : 1d-array

    The converted physical property values along the layer.

gradient(p)

The gradient vector of the misfit function.

\[\bar{g} = -2\bar{\bar{J}}^T\bar{r}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix and \(\bar{r}\) is the residual vector.

Parameters:

  • p
    : 1d-array

    The parameter vector where the gradient is evaluated

Returns:

  • gradient
    : 1d-array

    The gradient vector.

hessian(p)

The Hessian of the misfit function with respect to the parameters.

Calculated using the Gauss approximation:

\[\bar{\bar{H}} \approx 2\bar{\bar{J}}^T\bar{\bar{J}}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix.

For linear problems, the Hessian matrix is cached in memory, so calling this method again will not trigger a re-calculation.

Parameters:

  • p
    : 1d-array

    The parameter vector where the Hessian is evaluated

Returns:

  • hessian
    : 2d-array

    The Hessian matrix

jacobian(p)[source]

Calculate the Jacobian matrix for a given parameter vector.

predicted(p)

Calculate the data predicted by a given parameter vector.

Parameters:

  • p
    : 1d-array (optional)

    The parameter vector with the estimated physical properties of the layer. If not given, will use the value calculated by .fit().

Returns:

  • result
    : 1d-array

    The predicted data vector.

regul_param

The regularization parameter (scale factor) for the objetive function.

Defaults to 1.

residuals(p=None)

Calculate the residuals vector (observed - predicted data).

Parameters:

  • p
    : 1d-array or None

    The parameter vector used to calculate the residuals. If None, will use the current estimate stored in estimate_.

Returns:

  • residuals
    : 1d-array or list of 1d-arrays

    The residual vector. If this is the sum of 1 or more Misfit instances, will return the residual vector from each of the summed misfits in the order of the sum.

set_weights(weights)

Set the data weights.

Using weights for the data, the least-squares data-misfit function becomes:

\[\phi = \bar{r}^T \bar{\bar{W}}\bar{r}\]

Parameters:

  • weights
    : 1d-array or 2d-array or None

    Weights for the data vector. If None, will remove any weights that have been set before. If it is a 2d-array, it will be interpreted as the weight matrix \(\bar{\bar{W}}\). If it is a 1d-array, it will be interpreted as the diagonal of the weight matrix (all off-diagonal elements will default to zero). The weight matrix can be a sparse array from scipy.sparse.

value(p)

Calculate the value of the misfit for a given parameter vector.

The value is given by:

\[\phi(\bar{p}) = \bar{r}^T\bar{\bar{W}}\bar{r}\]

where \(\bar{r}\) is the residual vector and \(bar{\bar{W}}\) are optional data weights.

Parameters:

  • p
    : 1d-array or None

    The parameter vector.

Returns:

  • value
    : float

    The value of the misfit function.

class fatiando.gravmag.eqlayer.PELSmoothness(grid, windows, degree)[source]

Bases: fatiando.inversion.regularization.Smoothness

Regularization to “join” neighboring windows in the PEL.

Use this with PELGravity and PELTotalField.

Parameters passed to PELSmoothness must be the same as passed to the PEL solvers.

Parameters:

  • grid
    : PointGrid

    The sources in the equivalent layer.

  • windows
    : tuple = (ny, nx)

    The number of windows that the layer will be divided in the y and x directions, respectively.

  • degree
    : int

    The degree of the bivariate polynomials used in each window of the PEL

See the docstring of PELGravity for an example usage.

copy(deep=False)

Make a copy of me together with all the cached methods.

gradient(p)

Calculate the gradient vector.

Parameters:

  • p
    : 1d-array or None

    The parameter vector. If None, will return 0.

Returns:

  • gradient
    : 1d-array

    The gradient

hessian(p)

Calculate the Hessian matrix.

Parameters:

  • p
    : 1d-array

    The parameter vector

Returns:

  • hessian
    : 2d-array

    The Hessian

regul_param

The regularization parameter (scale factor) for the objetive function.

Defaults to 1.

value(p)

Calculate the value of this function.

Parameters:

  • p
    : 1d-array

    The parameter vector

Returns:

  • value
    : float

    The value of this function evaluated at p

class fatiando.gravmag.eqlayer.PELTotalField(x, y, z, data, inc, dec, grid, windows, degree, sinc=None, sdec=None)[source]

Bases: fatiando.gravmag.eqlayer.PELBase

Estimate a polynomial equivalent layer from magnetic total field anomaly.

Note

Assumes x = North, y = East, z = Down.

Parameters:

  • x, y, z
    : 1d-arrays

    The x, y, z coordinates of each data point.

  • data
    : 1d-array

    The total field magnetic anomaly data at each point.

  • inc, dec
    : floats

    The inclination and declination of the inducing field

  • grid
    : PointGrid

    The sources in the equivalent layer. Will invert for the magnetization intensity of each point in the grid.

  • windows
    : tuple = (ny, nx)

    The number of windows that the layer will be divided in the y and x directions, respectively

  • degree
    : int

    The degree of the bivariate polynomials used in each window of the PEL

  • sinc, sdec
    : None or floats

    The inclination and declination of the equivalent layer. Use these if there is remanent magnetization and the total magnetization of the layer if different from the induced magnetization. If there is only induced magnetization, use None

config(method, **kwargs)

Configure the optimization method and its parameters.

This sets the method used by fit and the keyword arguments that are passed to it.

Parameters:

  • method
    : string

    The optimization method. One of: 'linear', 'newton', 'levmarq', 'steepest', 'acor'

Other keyword arguments that can be passed are the ones allowed by each method.

Some methods have required arguments:

  • newton, levmarq and steepest require the initial argument (an initial estimate for the gradient descent)
  • acor requires the bounds argument (min/max values for the search space)

See the corresponding docstrings for more information:

copy(deep=False)

Make a copy of me together with all the cached methods.

estimate_

A nicely formatted version of the estimate.

If the class implements a fmt_estimate method, this will its results. This can be used to convert the parameter vector to a more useful form, like a fatiando.mesher object.

fit()

Solve for the parameter vector that minimizes this objective function.

Uses the optimization method and parameters defined using the config method.

The estimated parameter vector can be accessed through the p_ attribute. A (possibly) formatted version (converted to a more manageable type) of the estimate can be accessed through the property estimate_.

fmt_estimate(coefs)

Convert the estimated polynomial coefficients to physical property values along the layer.

Parameters:

  • coefs
    : 1d-array

    The estimated parameter vector with the polynomial coefficients

Returns:

  • estimate
    : 1d-array

    The converted physical property values along the layer.

gradient(p)

The gradient vector of the misfit function.

\[\bar{g} = -2\bar{\bar{J}}^T\bar{r}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix and \(\bar{r}\) is the residual vector.

Parameters:

  • p
    : 1d-array

    The parameter vector where the gradient is evaluated

Returns:

  • gradient
    : 1d-array

    The gradient vector.

hessian(p)

The Hessian of the misfit function with respect to the parameters.

Calculated using the Gauss approximation:

\[\bar{\bar{H}} \approx 2\bar{\bar{J}}^T\bar{\bar{J}}\]

where \(\bar{\bar{J}}\) is the Jacobian matrix.

For linear problems, the Hessian matrix is cached in memory, so calling this method again will not trigger a re-calculation.

Parameters:

  • p
    : 1d-array

    The parameter vector where the Hessian is evaluated

Returns:

  • hessian
    : 2d-array

    The Hessian matrix

jacobian(p)[source]

Calculate the Jacobian matrix for a given parameter vector.

predicted(p)

Calculate the data predicted by a given parameter vector.

Parameters:

  • p
    : 1d-array (optional)

    The parameter vector with the estimated physical properties of the layer. If not given, will use the value calculated by .fit().

Returns:

  • result
    : 1d-array

    The predicted data vector.

regul_param

The regularization parameter (scale factor) for the objetive function.

Defaults to 1.

residuals(p=None)

Calculate the residuals vector (observed - predicted data).

Parameters:

  • p
    : 1d-array or None

    The parameter vector used to calculate the residuals. If None, will use the current estimate stored in estimate_.

Returns:

  • residuals
    : 1d-array or list of 1d-arrays

    The residual vector. If this is the sum of 1 or more Misfit instances, will return the residual vector from each of the summed misfits in the order of the sum.

set_weights(weights)

Set the data weights.

Using weights for the data, the least-squares data-misfit function becomes:

\[\phi = \bar{r}^T \bar{\bar{W}}\bar{r}\]

Parameters:

  • weights
    : 1d-array or 2d-array or None

    Weights for the data vector. If None, will remove any weights that have been set before. If it is a 2d-array, it will be interpreted as the weight matrix \(\bar{\bar{W}}\). If it is a 1d-array, it will be interpreted as the diagonal of the weight matrix (all off-diagonal elements will default to zero). The weight matrix can be a sparse array from scipy.sparse.

value(p)

Calculate the value of the misfit for a given parameter vector.

The value is given by:

\[\phi(\bar{p}) = \bar{r}^T\bar{\bar{W}}\bar{r}\]

where \(\bar{r}\) is the residual vector and \(bar{\bar{W}}\) are optional data weights.

Parameters:

  • p
    : 1d-array or None

    The parameter vector.

Returns:

  • value
    : float

    The value of the misfit function.

fatiando.gravmag.eqlayer.ncoeffs(degree)[source]

Calculate the number of coefficients in a bivarite polynomail.

Parameters:

  • degree
    : int

    The degree of the polynomial

Returns:

  • n
    : int

    The number of coefficients

Examples:

>>> ncoeffs(1)
3
>>> ncoeffs(2)
6
>>> ncoeffs(3)
10
>>> ncoeffs(4)
15