geoist.pfm package

Subpackages

• geoist.pfm.cof package
  • Submodules
  • geoist.pfm.cof.noaamm module
  • Module contents
• geoist.pfm.test package
  • Submodules
  • geoist.pfm.test.test_grdio module
  • Module contents

Submodules

geoist.pfm.basin2d module

geoist.pfm.euler module

• EulerDeconv: The classic 3D solution to Euler's equation for potential fields (Reid et al., 1990). Runs on the whole dataset.

• EulerDeconvEW: Run Euler deconvolution on an expanding window over the data set and keep the best estimate.

• EulerDeconvMW: Run Euler deconvolution on a moving window over the data set to produce a set of estimates.

References

Reid, A. B., J. M. Allsop, H. Granser, A. J. Millett, and I. W. Somerton (1990), Magnetic interpretation in three dimensions using Euler deconvolution, Geophysics, 55(1), 80-91, doi:10.1190/1.1442774.

class geoist.pfm.euler.EulerDeconv(x, y, z, field, xderiv, yderiv, zderiv, structural_index)

基类：geoist.inversion.misfit.Misfit

Classic 3D Euler deconvolution of potential field data.

Follows the formulation of Reid et al. (1990). Performs the deconvolution on the whole data set. For windowed approaches, use EulerDeconvMW (moving window) and EulerDeconvEW (expanding window).

Works on any potential field that satisfies Euler's homogeneity equation (both gravity and magnetic, assuming simple sources):

\[(x_i - x_0)\dfrac{\partial f_i}{\partial x} + (y_i - y_0)\dfrac{\partial f_i}{\partial y} + (z_i - z_0)\dfrac{\partial f_i}{\partial z} = \eta (b - f_i),\]

in which \(f_i\) is the given potential field observation at point \((x_i, y_i, z_i)\), \(b\) is the base level (a constant shift of the field, like a regional field), \(\eta\) is the structural index, and \((x_0, y_0, z_0)\) are the coordinates of a point on the source (for a sphere, this is the center point).

The Euler deconvolution estimates \((x_0, y_0, z_0)\) and \(b\) given a potential field and its x, y, z derivatives and the structural index. However, this assumes that the sources are ideal (see the table below). We recommend reading Reid and Thurston (2014) for a discussion on what the structural index means and what it does not mean.

警告

Please read the paper Reid et al. (2014) to avoid doing horrible things with Euler deconvolution. Uieda et al. (2014) offer a practical tutorial using geoist code and show some common misinterpretations.

After Reid et al. (2014), values of the structural index (SI) can be:

Source type	SI (Mag)	SI (Grav)
Point, sphere	3	2
Line, cylinder, thin bed fault	2	1
Thin sheet edge, thin sill, thin dyke	1	0

Use the fit method to run the deconvolution. The estimated coordinates \((x_0, y_0, z_0)\) are stored in the estimate_ attribute and the estimated base level \(b\) is stored in baselevel_.

注解

Using structural index of 0 is not supported yet.

注解

The data does not need to be gridded for this! So long as you can calculate the derivatives of non-gridded data (using an Equivalent Layer, for example).

注解

x is North, y is East, and z is down.

注解

Units of the input data (x, y, z, field, derivatives) must be in SI units! Otherwise, the results will be in strange units. Use functions in geoist.pfm.giutils to convert between units.

Parameters:

• x, y, z1d-arrays

  The x, y, and z coordinates of the observation points

• field1d-array

  The potential field measured at the observation points

• xderiv, yderiv, zderiv1d-arrays

  The x-, y-, and z-derivatives of the potential field (measured or calculated) at the observation points

• indexfloat

  The structural index of the source

References:

Reid, A. B., J. M. Allsop, H. Granser, A. J. Millett, and I. W. Somerton (1990), Magnetic interpretation in three dimensions using Euler deconvolution, Geophysics, 55(1), 80-91, doi:10.1190/1.1442774.

Reid, A. B., J. Ebbing, and S. J. Webb (2014), Avoidable Euler Errors – the use and abuse of Euler deconvolution applied to potential fields, Geophysical Prospecting, doi:10.1111/1365-2478.12119.

Reid, A., and J. Thurston (2014), The structural index in gravity and magnetic interpretation: Errors, uses, and abuses, GEOPHYSICS, 79(4), J61-J66, doi:10.1190/geo2013-0235.1.

Uieda, L., V. C. Oliveira Jr., and V. C. F. Barbosa (2014), Geophysical tutorial: Euler deconvolution of potential-field data, The Leading Edge, 33(4), 448-450, doi:10.1190/tle33040448.1.

property baselevel_

fmt_estimate(p)

Separate the (x, y, z) point coordinates from the baselevel.

Coordinates are stored in estimate_ and a base level is stored in baselevel_.

jacobian(p)

predicted(p)

Calculate the predicted data for a given parameter vector.

Parameters:

• p1d-array or None

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

Returns:

• predicted1d-array or list of 1d-arrays

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

class geoist.pfm.euler.EulerDeconvEW(x, y, z, field, xderiv, yderiv, zderiv, structural_index, center, sizes)

基类：geoist.pfm.euler.EulerDeconv

Euler deconvolution using an expanding window scheme.

Uses data inside a window of growing size to perform the Euler deconvolution. Keeps the best result, judged by the estimated error.

The deconvolution is performed as in EulerDeconv.

Use the fit method to produce an estimate. The estimated point is stored in the attribute estimate_ and the base level in baselevel_.

Parameters:

• x, y, z1d-arrays

  The x, y, and z coordinates of the observation points

• field1d-array

  The potential field measured at the observation points

• xderiv, yderiv, zderiv1d-arrays

  The x-, y-, and z-derivatives of the potential field (measured or calculated) at the observation points

• indexfloat

  The structural index of the source

• center[x, y]

  The x, y coordinates of the center of the expanding windows.

• sizeslist or 1d-array

  The sizes of the windows.

fit()

Perform the Euler deconvolution with expanding windows.

The estimated point is stored in estimate_, the base level in baselevel_.

class geoist.pfm.euler.EulerDeconvMW(x, y, z, field, xderiv, yderiv, zderiv, structural_index, windows, size, keep=0.2)

基类：geoist.pfm.euler.EulerDeconv

Solve an Euler deconvolution problem using a moving window scheme.

Uses data inside a window moving to perform the Euler deconvolution. Keeps only a top percentage of the estimates from all windows.

The deconvolution is performed as in EulerDeconv.

Use the fit method to produce an estimate. The estimated points are stored in estimate_ as a 2D numpy array. Each line in the array is an [x, y, z] coordinate of a point. The base levels are stored in baselevel_.

Parameters:

• x, y, z1d-arrays

  The x, y, and z coordinates of the observation points

• field1d-array

  The potential field measured at the observation points

• xderiv, yderiv, zderiv1d-arrays

  The x-, y-, and z-derivatives of the potential field (measured or calculated) at the observation points

• indexfloat

  The structural index of the source

• windows(ny, nx)

  The number of windows in the y and x directions

• size(dy, dx)

  The size of the windows in the y and x directions

• keepfloat

  Decimal percentage of solutions to keep. Will rank the solutions by increasing error and keep only the first keep percent.

property baselevel_

fit()

Perform the Euler deconvolution on a moving window.

The estimated points are stored in estimate_, the base levels in baselevel_.

fmt_estimate(p)

Separate the (x, y, z) point coordinates from the baselevel.

Coordinates are stored in estimate_ and a base level is stored in baselevel_.

geoist.pfm.giconstants module

Holds all physical constants and unit conversions used in GEOIST. All modules should import the constants from here! All constants should be in SI, unless otherwise stated!

geoist.pfm.giconstants.CM = 1e-07

Proportionality constant used in the magnetic method in henry/m (SI)

geoist.pfm.giconstants.G = 6.673e-11

The gravitational constant in \(m^3 kg^{-1} s^{-1}\)

geoist.pfm.giconstants.MEAN_EARTH_RADIUS = 6378137.0

The mean earth radius in meters

geoist.pfm.giconstants.PERM_FREE_SPACE = 1.2566370614359173e-06

Permeability of free space in \(N A^{-2}\)

geoist.pfm.giconstants.SI2EOTVOS = 1000000000.0

\(1/s^2 = 10^9\ Eotvos\)

Type

Conversion factor from SI units to Eotvos

geoist.pfm.giconstants.SI2MGAL = 100000.0

\(1\ m/s^2 = 10^5\ mGal\)

Type

Conversion factor from SI units to mGal

geoist.pfm.giconstants.T2NT = 1000000000.0

Conversion factor from tesla to nanotesla

geoist.pfm.giconstants.THERMAL_DIFFUSIVITY = 1e-06

The default thermal diffusivity in \(m^2/s\)

geoist.pfm.giconstants.THERMAL_DIFFUSIVITY_YEAR = 31.5576

The default thermal diffusivity but in \(m^2/year\)

geoist.pfm.giutils module

Description : Miscellaneous utility functions.

class geoist.pfm.giutils.SparseList(size, elements=None)

基类：object

Store only non-zero elements on an immutable list.

Can iterate over and access elements just like if it were a list.

Parameters:

• sizeint

  Size of the list.

• elementsdict

  Dictionary used to initialize the list. Keys are the index of the elements and values are their respective values.

Example:

>>> l = SparseList(5)
>>> l[3] = 42.0
>>> print len(l)
5
>>> print l[1], l[3]
0.0 42.0
>>> l[1] += 3.0
>>> for i in l:
...     print i,
0.0 3.0 0.0 42.0 0.0
>>> l2 = SparseList(4, elements={1:3.2, 3:2.8})
>>> for i in l2:
...     print i,
0.0 3.2 0.0 2.8

next()

geoist.pfm.giutils.ang2vec(intensity, inc, dec)

Convert intensity, inclination and declination to a 3-component vector

注解

Coordinate system is assumed to be x->North, y->East, z->Down. Inclination is positive down and declination is measured with respect to x (North).

Parameter:

• intensityfloat or array

  The intensity (norm) of the vector

• incfloat

  The inclination of the vector (in degrees)

• decfloat

  The declination of the vector (in degrees)

Returns:

• vecarray = [x, y, z]

  The vector

Examples:

>>> import numpy
>>> print ang2vec(3, 45, 45)
[ 1.5         1.5         2.12132034]
>>> print ang2vec(numpy.arange(4), 45, 45)
[[ 0.          0.          0.        ]
 [ 0.5         0.5         0.70710678]
 [ 1.          1.          1.41421356]
 [ 1.5         1.5         2.12132034]]

geoist.pfm.giutils.check_hash(fname, known_hash, hash_type='sha256')

Assert that the hash of a file is equal to a known hash.

Useful for checking if a file has changed or been corrupted.

Parameters:

• fnamestring

  The name of the file.

• known_hashstring

  The known (recorded) hash of the file.

• hash_typestring

  What kind of hash is it. Can be anything defined in Python's hashlib.

Raises:

• AssertionError if the hash is different from the known hash.

geoist.pfm.giutils.contaminate(data, stddev, percent=False, return_stddev=False, seed=None)

Add pseudorandom gaussian noise to an array.

Noise added is normally distributed with zero mean.

Parameters:

• dataarray or list of arrays

  Data to contaminate

• stddevfloat or list of floats

  Standard deviation of the Gaussian noise that will be added to data

• percentTrue or False

  If True, will consider stddev as a decimal percentage and the standard deviation of the Gaussian noise will be this percentage of the maximum absolute value of data

• return_stddevTrue or False

  If True, will return also the standard deviation used to contaminate data

• seedNone or int

  Seed used to generate the pseudo-random numbers. If None, will use a different seed every time. Use the same seed to generate the same random sequence to contaminate the data.

Returns:

if return_stddev is False:

• contamarray or list of arrays

  The contaminated data array

else:

• resultslist = [contam, stddev]

  The contaminated data array and the standard deviation used to contaminate it.

Examples:

>>> import numpy as np
>>> data = np.ones(5)
>>> noisy = contaminate(data, 0.1, seed=0)
>>> print noisy
[ 1.03137726  0.89498775  0.95284582  1.07906135  1.04172782]
>>> noisy, std = contaminate(data, 0.05, seed=0, percent=True,
...                          return_stddev=True)
>>> print std
0.05
>>> print noisy
[ 1.01568863  0.94749387  0.97642291  1.03953067  1.02086391]
>>> data = [np.zeros(5), np.ones(3)]
>>> noisy = contaminate(data, [0.1, 0.2], seed=0)
>>> print noisy[0]
[ 0.03137726 -0.10501225 -0.04715418  0.07906135  0.04172782]
>>> print noisy[1]
[ 0.81644754  1.20192079  0.98163167]

geoist.pfm.giutils.dircos(inc, dec)

Returns the 3 coordinates of a unit vector given its inclination and declination.

注解

Coordinate system is assumed to be x->North, y->East, z->Down. Inclination is positive down and declination is measured with respect to x (North).

Parameter:

• incfloat

  The inclination of the vector (in degrees)

• decfloat

  The declination of the vector (in degrees)

Returns:

• vectlist = [x, y, z]

  The unit vector

geoist.pfm.giutils.eotvos2si(value)

Convert a value from Eotvos to SI units.

Parameters:

• valuenumber or array

  The value in Eotvos

Returns:

• valuenumber or array

  The value in SI

geoist.pfm.giutils.gaussian(x, mean, std)

Non-normalized Gaussian function

\[G(x,\bar{x},\sigma) = \exp\left(-\frac{(x-\bar{x})^2}{\sigma^2} \right)\]

Parameters:

• xfloat or array

  Values at which to calculate the Gaussian function

• meanfloat

  The mean of the distribution \(\bar{x}\)

• stdfloat

  The standard deviation of the distribution \(\sigma\)

Returns:

• gaussarray

  Gaussian function evaluated at x

geoist.pfm.giutils.gaussian2d(x, y, sigma_x, sigma_y, x0=0, y0=0, angle=0.0)

Non-normalized 2D Gaussian function

Parameters:

• x, yfloat or arrays

  Coordinates at which to calculate the Gaussian function

• sigma_x, sigma_yfloat

  Standard deviation in the x and y directions

• x0, y0float

  Coordinates of the center of the distribution

• anglefloat

  Rotation angle of the gaussian measure from the x axis (north) growing positive to the east (positive y axis)

Returns:

• gaussarray

  Gaussian function evaluated at x, y

geoist.pfm.giutils.mgal2si(value)

Convert a value from mGal to SI units.

Parameters:

• valuenumber or array

  The value in mGal

Returns:

• valuenumber or array

  The value in SI

geoist.pfm.giutils.nt2si(value)

Convert a value from nanoTesla to SI units.

Parameters:

• valuenumber or array

  The value in nanoTesla

Returns:

• valuenumber or array

  The value in SI

geoist.pfm.giutils.safe_diagonal(matrix)

Get the diagonal of a matrix using the appropriate method.

Parameters:

• matrix2d-array, matrix, sparse matrix

  The matrix...

Returns:

• diag1d-array

  A numpy array with the diagonal of the matrix

geoist.pfm.giutils.safe_dot(a, b)

Make the dot product using the appropriate method.

If a and b are dense, will use numpy.dot. If either is sparse (from scipy.sparse) will use the multiplication operator (i.e., *).

Parameters:

• a, barray or matrix

  The vectors/matrices to take the dot product of.

Returns:

• prodarray or matrix

  The dot product of a and b

geoist.pfm.giutils.safe_inverse(matrix)

Calculate the inverse of a matrix using an apropriate algorithm.

Uses the standard numpy.linalg.inv if matrix is dense. If it is sparse (from scipy.sparse) then will use scipy.sparse.linalg.inv.

Parameters:

• matrix2d-array

  The matrix

Returns:

• inverse2d-array

  The inverse of matrix

geoist.pfm.giutils.safe_solve(matrix, vector)

Solve a linear system using an apropriate algorithm.

Uses the standard numpy.linalg.solve if both matrix and vector are dense.

If any of the two is sparse (from scipy.sparse) then will use the Conjugate Gradient Method (scipy.sparse.cgs).

Parameters:

• matrix2d-array

  The matrix defining the linear system

• vector1d or 2d-array

  The right-side vector of the system

Returns:

• solution1d or 2d-array

  The solution of the linear system

geoist.pfm.giutils.si2eotvos(value)

Convert a value from SI units to Eotvos.

Parameters:

• valuenumber or array

  The value in SI

Returns:

• valuenumber or array

  The value in Eotvos

geoist.pfm.giutils.si2mgal(value)

Convert a value from SI units to mGal.

Parameters:

• valuenumber or array

  The value in SI

Returns:

• valuenumber or array

  The value in mGal

geoist.pfm.giutils.si2nt(value)

Convert a value from SI units to nanoTesla.

Parameters:

• valuenumber or array

  The value in SI

Returns:

• valuenumber or array

  The value in nanoTesla

geoist.pfm.giutils.sph2cart(lon, lat, height)

Convert spherical coordinates to Cartesian geocentric coordinates.

Parameters:

• lon, lat, heightfloats

  Spherical coordinates. lon and lat in degrees, height in meters. height is the height above mean Earth radius.

Returns:

• x, y, zfloats

  Converted Cartesian coordinates

geoist.pfm.giutils.vec2ang(vector)

Convert a 3-component vector to intensity, inclination and declination.

注解

Coordinate system is assumed to be x->North, y->East, z->Down. Inclination is positive down and declination is measured with respect to x (North).

Parameter:

• vectorarray = [x, y, z]

  The vector

Returns:

• [intensity, inclination, declination]floats

  The intensity, inclination and declination (in degrees)

Examples:

>>> s = vec2ang([1.5, 1.5, 2.121320343559643])
>>> print "%.3f %.3f %.3f" % tuple(s)
3.000 45.000 45.000

geoist.pfm.grdio module

Name : grdio.py Created on : 2018/11/24 08:57 Author : Steve Chen <chenshi@cea-igp.ac.cn> Affiliation : Institute of Geophysics, CEA. Version : 0.1.0 Copyright : Copyright (C) 2018-2020 GEOIST Development Team. All Rights Reserved. License : Distributed under the MIT License. See LICENSE.txt for more info. Github : https://igp-gravity.github.io/ Description : Application for processing grid data of potential dataset.

class geoist.pfm.grdio.grddata

基类：object

Grid Data Object .. attribute:: data

array to contain raster data

type

numpy masked array

xmin

min value X coordinate of raster grid

Type

float

ymin

min value Y coordinate of raster grid

Type

float

xdim

x-dimension of grid cell

Type

float

ydim

y-dimension of grid cell

Type

float

typeofdata

number of datatype

Type

int

dataname

data name or id

Type

str

rows

number of rows for each raster grid/band

Type

int

cols

number of columns for each raster grid/band

Type

int

nullvalue

grid null or nodata value

Type

float

norm

normalized data

Type

dictionary

gtr

projection information

Type

tuple

wkt

projection information

Type

str

units

description of units to be used with color bars

Type

str

export_ascii(fname)

Export Ascii file

参数

data (grid Data) -- dataset to export

export_surfer(fname, flag=True, file_format='binary')

Export a surfer grid

参数

• fname (filename of grid dataset to export) --

• flag (True - Output Grid Grid) -- False - Output Bak Grid Grid

• file_format (binary/b - output binary format) -- ascii/a - output ascii format

fill_nulls(method='nearest')

Fill in the NaNs or masked values on interpolated points using nearest neighbors. method='nearest' or 'linear' or 'cubic'

grd2xyz(flag=True)

Return x,y,z 1-D array data from 2-D grid array.

参数

flag -- True - Output Grid Grid False - Output Bak Grid Grid

返回

x,y,z 1-D array data

load_ascii(fname, dtype='float64')

Load Ascii file

参数

data (grid Data) -- dataset to export

load_grd(fname, *args, **kwargs)

load_surfer(fname, *args, **kwargs)

Read data from a Surfer grid file.

Parameters:

• fnamestr

  Name of the Surfer grid file

• dtypenumpy dtype object or string

  The type of variable used for the data. Default is numpy.float64 for ascii data and is '=f' for binary data. Use numpy.float32 if the data are large and precision is not an issue.

• header_formatheader format (excluding the leading 'DSBB') following

  the convention of the struct module. Only used for binary data.

Returns:

geoist.pfm.grdio.regular(area, shape, z=None)

Create a regular grid.

The x directions is North-South and y East-West. Imagine the grid as a matrix with x varying in the lines and y in columns.

Returned arrays will be flattened to 1D with numpy.ravel.

Parameters:

• area

  (x1, x2, y1, y2): Borders of the grid

• shape

  Shape of the regular grid, ie (nx, ny).

• z

  Optional. z coordinate of the grid points. If given, will return an array with the value z.

Returns:

• [x, y]

  Numpy arrays with the x and y coordinates of the grid points

• [x, y, z]

  If z given. Numpy arrays with the x, y, and z coordinates of the grid points

Examples:

>>> x, y = regular((0, 10, 0, 5), (5, 3))
>>> print(x)
[  0.    0.    0.    2.5   2.5   2.5   5.    5.    5.    7.5   7.5   7.5
  10.   10.   10. ]
>>> print(x.reshape((5, 3)))
[[  0.    0.    0. ]
 [  2.5   2.5   2.5]
 [  5.    5.    5. ]
 [  7.5   7.5   7.5]
 [ 10.   10.   10. ]]

geoist.pfm.grdio.spacing(area, shape)

Returns the spacing between grid nodes

Parameters:

• area

  (x1, x2, y1, y2): Borders of the grid

• shape

  Shape of the regular grid, ie (nx, ny).

Returns:

• [dx, dy]

  Spacing the y and x directions

Examples:

>>> print(spacing((0, 10, 0, 20), (11, 11)))
[1.0, 2.0]
>>> print(spacing((0, 10, 0, 20), (11, 21)))
[1.0, 1.0]
>>> print(spacing((0, 10, 0, 20), (5, 21)))
[2.5, 1.0]
>>> print(spacing((0, 10, 0, 20), (21, 21)))
[0.5, 1.0]

geoist.pfm.hawaii_gravity module

Load gravity data from the eigen-6c4 model for Hawaii.

geoist.pfm.hawaii_gravity.fetch_hawaii_gravity()

Load gravity and topography data for Hawaii.

The raw gravity and topography data were generated from the eigen-6c4 spherical harmonic model and downloaded from ICGEM).

The topography is in meters and gravity in mGal.

x (north-south) and y (east-west) coordinates are UTM zone 4 (WGS84) in meters.

The gravity disturbance was calculated using the closed-form formula and the WGS84 ellipsoid.

The topography-free disturbances (both Bouguer and full) used densities 2670 kg/m³ for the crust and 1040 kg/m³ for the ocean water.

The full topography-free disturbance was calculated using a tesseroid model of the topography.

Returns:

• datadict

  A dictionary with the data fields:

  • 'metadata'string

    The headers from the original ICGEM files.

  • 'shape'tuple

    (nlat, nlon), the number of points in each dimension of the data grid.

  • 'area'tuple

    (south, north, west, east), the bounding dimensions of the data grid in degrees.

  • 'height'float

    The observation height of the gravity data in meters.

  • 'lat', 'lon'1d arrays

    The latitude and longitude coordinates of each point.

  • 'x', 'y'1d arrays

    The UTM zone 4 coordinates of each point in meters. x is north-south and y is east-west.

  • 'topography'1d array

    The topographic (geometric) height of each point in meters.

  • 'gravity'1d array

    The raw gravity value of each point in mGal.

  • 'disturbance'1d array

    The gravity disturbance value of each point in mGal.

  • 'topo-free'1d array

    The topography-free gravity disturbance value of each point in mGal calculated using tesseroids.

  • 'topo-free-bouguer'1d array

    The topography-free gravity disturbance value of each point in mGal calculated using the Bouguer plate approximation.

geoist.pfm.igrf module

This program, originally written in FORTRAN, was developed using subroutines

written by : A. Zunde

USGS, MS 964, Box 25046 Federal Center, Denver, Co. 80225

and

S.R.C. Malin & D.R. Barraclough

Institute of Geological Sciences, United Kingdom.

Translated

into C by : Craig H. Shaffer

29 July, 1988

Rewritten by : David Owens

For Susan McLean

Maintained by: Stefan Maus

Contact : stefan.maus@noaa.gov

National Geophysical Data Center

World Data Center-A for Solid Earth Geophysics

NOAA, E/GC1, 325 Broadway,

Boulder, CO 80303

class geoist.pfm.igrf.IGRF(parent=None)

基类：object

IGRF field calculation

This produces two datasets. The first is an IGRF dataset for the area of interest, defined by some input magnetic dataset. The second is the IGRF corrected form of that input magnetic dataset.

To do this, the input dataset must be reprojected from its local projection to degrees, where the IGRF correction will take place. This is done within this class.

parent

reference to the parent routine

Type

parent

indata

dictionary of input datasets

Type

dictionary

outdata

dictionary of output datasets

Type

dictionary

参数

• altmin (Double) -- Minimum height of selected model.

• altmax (Double array) -- array of MAXMOD Maximum height of model.

• maxalt (Double) -- Maximum height of selected model.

• d (float) -- Declination of the field from the geographic north (deg).

• sdate (float) -- start date inputted

• ddot (float) -- annual rate of change of decl. (arc-min/yr)

• alt (float) -- altitude above WGS84 Ellipsoid

• epoch (list) -- list of MAXMOD epoch of model.

• latitude (float) -- Latitude.

• longitude (float) -- Longitude.

• gh (numpy array) -- Schmidt quasi-normal internal spherical harmonic coeff. Schmidt quasi-normal internal spherical harmonic coeff. Coefficients of resulting model. Coefficients of rate of change model.

• i (float) -- Inclination (deg).

• idot (float) -- Rate of change of i (arc-min/yr).

• igdgc (int) -- Flag for geodetic or geocentric coordinate choice.

• irec_pos (int array) -- array of MAXMOD Record counter for header

• fileline (int) -- Current line in file (for errors)

• max1 (list, int) -- array of MAXMOD Main field coefficient.

• max2 (list, int) -- array of MAXMOD Secular variation coefficient.

• max3 (list, int) -- array of MAXMOD Acceleration coefficient.

• minyr (float) -- Min year of all models

• maxyr (float) -- Max year of all models

• yrmax (list, float) -- array of MAXMOD Max year of model.

• yrmin (list, float) -- array of MAXMOD Min year of model.

dihf(gh)

Computes the geomagnetic d, i, h, and f from x, y, and z.

FORTRAN : A. Zunde, USGS, MS 964, box 25046 Federal Center, Denver,

CO. 80225

C : C. H. Shaffer, Lockheed Missiles and Space Company, Sunnyvale CA

参数

• x (float) -- northward component

• y (float) -- eastward component

• z (float) -- vertically-downward component

返回

• d (float) -- declination

• i (float) -- inclination

• h (float) -- horizontal intensity

• f (float) -- total intensity

extrapsh(date, dte1, nmax1, nmax2, gh)

Extrapolates linearly a spherical harmonic model with a rate-of-change model.

FORTRAN : A. Zunde, USGS, MS 964, box 25046 Federal Center, Denver,

CO. 80225

C : C. H. Shaffer, Lockheed Missiles and Space Company, Sunnyvale CA

参数

• date (float) -- date of resulting model (in decimal year)

• dte1 (float) -- date of base model

• nmax1 (int) -- maximum degree and order of base model

• gh (numpy array) -- Schmidt quasi-normal internal spherical harmonic coefficients of base model and rate-of-change model

• nmax2 (int) -- maximum degree and order of rate-of-change model

返回

• gh (numpy array) -- Schmidt quasi-normal internal spherical harmonic coefficients

• nmax (int) -- maximum degree and order of resulting model

getshc(file, iflag, strec, nmax_of_gh, gh)

Reads spherical harmonic coefficients from the specified model into an array.

FORTRAN: Bill Flanagan, NOAA CORPS, DESDIS, NGDC, 325 Broadway,

Boulder CO. 80301

C: C. H. Shaffer, Lockheed Missiles and Space Company, Sunnyvale CA

参数

• file (file) -- reference to a file object

• iflag -- Flag for SV equal to 1 or not equal to 1 for designated read statements

• strec (int) -- Starting record number to read from model

• nmax_of_gh (int) -- Maximum degree and order of model

返回

gh -- Schmidt quasi-normal internal spherical harmonic coefficients

返回类型

list

interpsh(date, dte1, nmax1, dte2, nmax2, gh)

Interpolates linearly, in time, between two spherical harmonic models.

FORTRAN : A. Zunde, USGS, MS 964, box 25046 Federal Center, Denver,

CO. 80225

C : C. H. Shaffer, Lockheed Missiles and Space Company, Sunnyvale CA

参数

• date (float) -- date of resulting model (in decimal year)

• dte1 (float) -- date of earlier model

• nmax1 (int) -- maximum degree and order of earlier model

• gh (numpy array) -- Schmidt quasi-normal internal spherical harmonic coefficients of earlier model and internal model

• dte2 (float) -- date of later model

• nmax2 (int) -- maximum degree and order of later model

返回

• gh (numpy array) -- coefficients of resulting model

• nmax (int) -- maximum degree and order of resulting model

pnt(latitude, longitude, alt, magdata, Model='IGRF12')

Settings Dialog. This is the main entrypoint into this routine. It also contains the main IGRF code.

shval3(igdgc, flat, flon, elev, nmax, gh)

Calculates field components from spherical harmonic (sh) models.

Based on subroutine 'igrf' by D. R. Barraclough and S. R. C. Malin, report no. 71/1, institute of geological sciences, U.K.

FORTRAN : Norman W. Peddie, USGS, MS 964, box 25046 Federal Center,

Denver, CO. 80225

C : C. H. Shaffer, Lockheed Missiles and Space Company, Sunnyvale CA

参数

• igdgc (int) -- indicates coordinate system used set equal to 1 if geodetic, 2 if geocentric

• latitude (float) -- north latitude, in degrees

• longitude (float) -- east longitude, in degrees

• elev (float) -- WGS84 altitude above ellipsoid (igdgc=1), or radial distance from earth's center (igdgc=2)

• a2,b2 (float) -- squares of semi-major and semi-minor axes of the reference spheroid used for transforming between geodetic and geocentric coordinates or components

• nmax (int) -- maximum degree and order of coefficients

返回

• x (float) -- northward component

• y (float) -- eastward component

• z (float) -- vertically-downward component

geoist.pfm.imaging module

geoist.pfm.inv3d module

Funtion:3D grave&mag inversion using the space-domain approach with Bayes inference

CLASS list:

Model ModelBuilder ModelTs ModelTsCollection Source HexSource Field GravField MagField Mesh RegularMesh IrregulerMesh

Using packages:

~geoist.inversion.abic.py ~geoist.inversion.walsh.py ~geoist.inversion.toeplitz.py

class geoist.pfm.inv3d.Dimension(nx, ny, nz)

基类：tuple

property nx

Alias for field number 0

property ny

Alias for field number 1

property nz

Alias for field number 2

geoist.pfm.inv3d.Dimension2

geoist.pfm.inv3d.Dimension 的别名

class geoist.pfm.inv3d.Field

基类：object

ftype = ''

get_mesh()

get_name()

get_type()

get_unit()

mesh = ''

name = ''

unit = ''

class geoist.pfm.inv3d.GravField(name, unit, mesh)

基类：geoist.pfm.inv3d.Field

gen_data()

class geoist.pfm.inv3d.HexSource(name, unit='km')

基类：geoist.pfm.inv3d.Source

dims = None

gen_single_body(props=1.0)

load_from_txt()

props = None

show()

vols = None

class geoist.pfm.inv3d.IrregulerMesh(name, x=[], y=[], values=[])

基类：geoist.pfm.inv3d.Mesh

load_mesh_from_txt()

show()

values = []

x = []

y = []

class geoist.pfm.inv3d.MagField(name, unit, mesh)

基类：geoist.pfm.inv3d.Field

class geoist.pfm.inv3d.Margin(left, right, down, up)

基类：tuple

property down

Alias for field number 2

property left

Alias for field number 0

property right

Alias for field number 1

property up

Alias for field number 3

class geoist.pfm.inv3d.Mesh

基类：object

name = ''

class geoist.pfm.inv3d.Model(ModelBuilder, name)

基类：object

the interface class for inversion with InvModel in pfmodel.py

add_ref_stack(ref_model)

depth_constraint = None

field = ''

forward()

gen_kernel()

inverse()

load_from_json(filename)

load model parameters and initialize it from JSON file TO_DO LIST 1:BEI

margin = None

name = ''

obs_stack = []

path = ''

ref_stack = []

run(mean=False)

save_to_json(filename)

Export instance of model class to JSON file

set_datapath(data_path)

set_depth_constraint(kind=1)

set_margin()

set_refmod()

set_smooth(smooth)

set_weight(weight)

smooth = None

source = ''

weight = None

class geoist.pfm.inv3d.ModelBuilder(use_gpu=0)

基类：object

add_obs(obs)

add_ref(ref)

bfield = ''

bsource = ''

builder(name='demo')

logging = []

set_field(field)

set_problem(cell='prism', problem='3D')

set_source(source, dims, vols)

class geoist.pfm.inv3d.ModelTs(ModelBuilder, name)

基类：geoist.pfm.inv3d.Model

optimize_weights = ['obs', 'dx', 'dt']

run()

save_to_json(filename)

Export instance of model class to JSON file

set_smooth(smooth_components)

weights = {'dt': 200, 'dx': 1, 'dxy': 1, 'obs': 10}

class geoist.pfm.inv3d.ModelsCollection(ModelBuilder, name)

基类：geoist.pfm.inv3d.ModelTs

class geoist.pfm.inv3d.RegularMesh(name, xset=[], yset=[], values=[])

基类：geoist.pfm.inv3d.Mesh

load_mesh_from_grd()

show()

values = []

xset = []

yset = []

class geoist.pfm.inv3d.Smooth(dx, dy, dxy, dz, dt)

基类：tuple

property dt

Alias for field number 4

property dx

Alias for field number 0

property dxy

Alias for field number 2

property dy

Alias for field number 1

property dz

Alias for field number 3

class geoist.pfm.inv3d.Source

基类：object

ftype = ''

get_name()

get_type()

get_unit()

name = ''

unit = ''

class geoist.pfm.inv3d.Volume(xmin, xmax, ymin, ymax, zmin, zmax)

基类：tuple

property xmax

Alias for field number 1

property xmin

Alias for field number 0

property ymax

Alias for field number 3

property ymin

Alias for field number 2

property zmax

Alias for field number 5

property zmin

Alias for field number 4

class geoist.pfm.inv3d.Weight(dx, dy, dz, depth, obs, refer, bound)

基类：tuple

property bound

Alias for field number 6

property depth

Alias for field number 3

property dx

Alias for field number 0

property dy

Alias for field number 1

property dz

Alias for field number 2

property obs

Alias for field number 4

property refer

Alias for field number 5

geoist.pfm.magdir module

Estimation of the total magnetization vector of homogeneous bodies.

It estimates parameters related to the magnetization vector of homogeneous bodies.

Algorithms

• DipoleMagDir: This class estimates the Cartesian components of the magnetization vector of homogeneous dipolar bodies with known center. The estimated magnetization vector is converted to dipole moment, inclination (positive down) and declination (with respect to x, North).

---

class geoist.pfm.magdir.DipoleMagDir(x, y, z, data, inc, dec, points)

基类：geoist.inversion.misfit.Misfit

Estimate the magnetization vector of a set of dipoles from magnetic total field anomaly using the method of Oliveira Jr. et al. (2015).

By using the well-known first-order approximation of the total field anomaly (Blakely, 1996, p. 179) produced by a set of dipoles, the estimation of the Cartesian components of the magnetization vectors is formulated as linear inverse problem. After estimating the magnetization vectors, they are converted to dipole moment, inclination (positive down) and declination (with respect to x, North).

After solving, use the estimate_ attribute to get the estimated magnetization vectors in dipole moment, inclination and declination. The estimated magnetization vectors in Cartesian coordinates can be accessed through the p_ attribute.

注解

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

Parameters:

• x, y, z1d-arrays

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

• data1d-array

  The total field magnetic anomaly data at each point.

• inc, decfloats

  The inclination and declination of the inducing field

• pointslist of points [x, y, z]

  Each point [x, y, z] is the center of a dipole. Will invert for the Cartesian components of the magnetization vector of each dipole. Subsequently, the estimated magnetization vectors are converted to dipole moment, inclination and declination.

注解

Inclination is positive down and declination is measured with respect to x (North).

References:

Blakely, R. (1996), Potential theory in gravity and magnetic applications: CUP

Oliveira Jr., V. C., D. P. Sales, V. C. F. Barbosa, and L. Uieda (2015), Estimation of the total magnetization direction of approximately spherical bodies, Nonlin. Processes Geophys., 22(2), 215-232, doi:10.5194/npg-22-215-2015.

Examples:

Estimation of the total magnetization vector of dipoles with known centers

>>> import numpy as np
>>> from fatiando import gridder, utils
>>> from fatiando.gravmag import sphere
>>> from fatiando.mesher import Sphere, Prism
>>> # Produce some synthetic data
>>> area = (0, 10000, 0, 10000)
>>> x, y, z = gridder.scatter(area, 500, z=-150, seed=0)
>>> model = [Sphere(3000, 3000, 1000, 1000,
...              {'magnetization': utils.ang2vec(6.0, -20.0, -10.0)}),
...          Sphere(7000, 7000, 1000, 1000,
...              {'magnetization': utils.ang2vec(6.0, 30.0, -40.0)})]
>>> inc, dec = -9.5, -13
>>> tf = sphere.tf(x, y, z, model, inc, dec)
>>> # Give the coordinates of the dipoles
>>> points = [[3000.0, 3000.0, 1000.0], [7000.0, 7000.0, 1000.0]]
>>> p_true = np.hstack((ang2vec(CM*(4.*np.pi/3.)*6.0*1000**3,
...                                             -20.0, -10.0),
...                        ang2vec(CM*(4.*np.pi/3.)*6.0*1000**3,
...                                              30.0, -40.0)))
>>> estimate_true = [utils.vec2ang(p_true[3*i : 3*i + 3]) for i
...                                in range(len(points))]
>>> # Make a solver and fit it to the data
>>> solver = DipoleMagDir(x, y, z, tf, inc, dec, points).fit()
>>> # Check the fit
>>> np.allclose(tf, solver.predicted(), rtol=0.001, atol=0.001)
True
>>> # solver.p_ returns the Cartesian components of the
>>> # estimated magnetization vectors
>>> for p in solver.p_: print "%.10f" % p
2325.8255393651
-410.1057950109
-859.5903757213
1667.3411086852
-1399.0653093445
1256.6370614359
>>> # Check the estimated parameter vector
>>> np.allclose(p_true, solver.p_, rtol=0.001, atol=0.001)
True
>>> # The parameter vector is not that useful so use solver.estimate_
>>> # to convert the estimated magnetization vectors in dipole moment,
>>> # inclination and declination.
>>> for e in solver.estimate_:
...    print "%.10f %.10f %.10f" % (e[0], e[1], e[2])
2513.2741228718 -20.0000000000 -10.0000000000
2513.2741228718 30.0000000000 -40.0000000000
>>> # Check the converted estimate
>>> np.allclose(estimate_true, solver.estimate_, rtol=0.001,
...                                                 atol=0.001)
True

fmt_estimate(p)

Convert the estimate parameters from Cartesian to inclination, declication, and intensity.

jacobian(p)

predicted(p)

Calculate the predicted data for a given parameter vector.

Parameters:

• p1d-array or None

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

Returns:

• predicted1d-array or list of 1d-arrays

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

geoist.pfm.normgra module

Reference ellipsoids

This module uses instances of ReferenceEllipsoid to store the physical constants of ellipsoids. To create a new ellipsoid, just instantiate ReferenceEllipsoid and give it the semimajor axis a, the flattening f, the geocentric gravitational constant GM, and the angular velocity omega. All other quantities, like the gravity at the poles and equator, eccentricities, etc, are computed by the class from these 4 parameters.

Available ellipsoids:

• WGS84 (values taken from Hofmann-Wellenhof and Moritz, 2005):

  >>> from geoist.pfm.normgra import WGS84
  >>> print(WGS84.name)
  World Geodetic System 1984
  >>> print('{:.0f}'.format(WGS84.a))
  6378137
  >>> print('{:.17f}'.format(WGS84.f))
  0.00335281066474748
  >>> print('{:.10g}'.format(WGS84.GM))
  3.986004418e+14
  >>> print('{:.7g}'.format(WGS84.omega))
  7.292115e-05
  >>> print('{:.4f}'.format(WGS84.b))
  6356752.3142
  >>> print('{:.8f}'.format(WGS84.E)) # Linear eccentricity
  521854.00842339
  >>> print('{:.15f}'.format(WGS84.e_prime)) # second eccentricity
  0.082094437949696
  >>> print('{:.10f}'.format(WGS84.gamma_a)) # gravity at the equator
  9.7803253359
  >>> print('{:.11f}'.format(WGS84.gamma_b)) # gravity at the pole
  9.83218493786
  >>> print('{:.14f}'.format(WGS84.m))
  0.00344978650684
  

Normal gravity

• gamma_somigliana: compute the normal gravity using the Somigliana formula (Hofmann-Wellenhof and Moritz, 2005). Calculated on the ellipsoid.

• gamma_somigliana_free_air: compute normal gravity at a height using the Somigliana formula plus the free-air correction \(-0.3086H\ mGal/m\).

• gamma_closed_form: compute normal gravity using the closed form expression from Li and Gotze (2001). Can compute anywhere (on, above, under the ellipsoid).

Bouguer

• bouguer_plate: compute the gravitational attraction of an infinite plate (Bouguer plate). Calculated on top of the plate.

You can use pfm.prism and pfm.tesseroid to calculate the terrain effect for a better correction.

References

Hofmann-Wellenhof, B. and H. Moritz, 2005, Physical Geodesy, Springer-Verlag Wien, ISBN-13: 978-3-211-23584-3

Li, X. and H. J. Gotze, 2001, Tutorial: Ellipsoid, geoid, gravity, geodesy, and geophysics, Geophysics, 66(6), p. 1660-1668, doi: 10.1190/1.1487109

---

class geoist.pfm.normgra.ReferenceEllipsoid(name, a, f, GM, omega)

基类：object

A generic reference ellipsoid.

It stores the physical constants defining the ellipsoid and has properties for computing other (derived) quantities.

All quantities are expected and returned in SI units.

Parameters:

• afloat

  The semimajor axis (the largest one, at the equator). In meters.

• ffloat

  The flattening. Adimensional.

• GMfloat

  The geocentric gravitational constant of the earth, including the atmosphere. In \(m^3 s^{-2}\).

• omegafloat

  The angular velocity of the earth. In \(rad s^{-1}\).

property E

Linear eccentricity

property GM

Geocentric gravitational constant (including the atmosphere)

property a

Semimajor axis

property b

Semiminor axis

property e_prime

Second eccentricity

property f

Flattening

property gamma_a

Normal gravity at the equator

property gamma_b

Normal gravity at the poles

property m

\(\omega^2 a^2 b / (GM)\)

property omega

Angular velocity

geoist.pfm.normgra.bouguer_plate(topography, density_rock=2670, density_water=1040)

Calculate the gravitational effect of an infinite Bouguer plate.

注解

The effect is calculated on top of the plate.

Uses the famous \(g_{BG} = 2 \pi G \rho t\) formula, where t is the height of the topography. On water (i.e., t < 0), uses \(g_{BG} = 2 \pi G (\rho_{water} - \rho_{rock})\times (-t)\).

Parameters:

• topographyfloat or numpy array

  The height of topography (in meters).

• density_rockfloat

  The density of crustal rocks

• density_waterfloat

  The density of ocean water

Returns:

• g_bouguerfloat or array

  The computed gravitational effect of the Bouguer plate

geoist.pfm.normgra.gamma_closed_form(latitude, height, ellipsoid=<geoist.pfm.normgra.ReferenceEllipsoid object>)

Calculate normal gravity at a height using the closed form expression of Li and Gotze (2001).

Parameters:

• latitudefloat or numpy array

  The latitude where the normal gravity will be computed (in degrees)

• heightfloat or numpy array

  The height of computation (in meters). Should be ellipsoidal (geometric) heights for geophysical purposes.

• ellipsoidReferenceEllipsoid

  The reference ellipsoid used.

Returns:

• gammafloat or numpy array

  The computed normal gravity (in mGal).

geoist.pfm.normgra.gamma_somigliana(latitude, ellipsoid=<geoist.pfm.normgra.ReferenceEllipsoid object>)

Calculate the normal gravity by using Somigliana's formula.

This formula computes normal gravity on the ellipsoid (height = 0).

Parameters:

• latitudefloat or numpy array

  The latitude where the normal gravity will be computed (in degrees)

• ellipsoidReferenceEllipsoid

  The reference ellipsoid used.

Returns:

• gammafloat or numpy array

  The computed normal gravity (in mGal).

geoist.pfm.normgra.gamma_somigliana_free_air(latitude, height, ellipsoid=<geoist.pfm.normgra.ReferenceEllipsoid object>)

Calculate the normal gravity at a height using Somigliana's formula and the free-air correction.

Parameters:

• latitudefloat or numpy array

  The latitude where the normal gravity will be computed (in degrees)

• heightfloat or numpy array

  The height of computation (in meters). Should be ellipsoidal (geometric) heights for geophysical purposes.

• ellipsoidReferenceEllipsoid

  The reference ellipsoid used.

Returns:

• gammafloat or numpy array

  The computed normal gravity (in mGal).

geoist.pfm.pftrans module

Potential field transformations, like upward continuation and derivatives. .. note:: Most, if not all, functions here required gridded data.

geoist.pfm.pftrans.derivx(x, y, data, shape, order=1, method='fd')

Calculate the derivative of a potential field in the x direction.

注解

Requires gridded data.

警告

If the data is not in SI units, the derivative will be in strange units! I strongly recommend converting the data to SI before calculating the derivative (use one of the unit conversion functions of geoist.utils). This way the derivative will be in SI units and can be easily converted to what unit you want.

Parameters:

• x, y1D-arrays

  The x and y coordinates of the grid points

• data1D-array

  The potential field at the grid points

• shapetuple = (nx, ny)

  The shape of the grid

• orderint

  The order of the derivative

• methodstring

  The method used to calculate the derivatives. Options are: 'fd' for central finite-differences (more stable) or 'fft' for the Fast Fourier Transform.

Returns:

• deriv1D-array

  The derivative

geoist.pfm.pftrans.derivy(x, y, data, shape, order=1, method='fd')

Calculate the derivative of a potential field in the y direction.

注解

Requires gridded data.

警告

If the data is not in SI units, the derivative will be in strange units! I strongly recommend converting the data to SI before calculating the derivative (use one of the unit conversion functions of geoist.utils). This way the derivative will be in SI units and can be easily converted to what unit you want.

Parameters:

• x, y1D-arrays

  The x and y coordinates of the grid points

• data1D-array

  The potential field at the grid points

• shapetuple = (nx, ny)

  The shape of the grid

• orderint

  The order of the derivative

• methodstring

  The method used to calculate the derivatives. Options are: 'fd' for central finite-differences (more stable) or 'fft' for the Fast Fourier Transform.

Returns:

• deriv1D-array

  The derivative

geoist.pfm.pftrans.derivz(x, y, data, shape, order=1, method='fft')

Calculate the derivative of a potential field in the z direction.

注解

Requires gridded data.

警告

If the data is not in SI units, the derivative will be in strange units! I strongly recommend converting the data to SI before calculating the derivative (use one of the unit conversion functions of geoist.utils). This way the derivative will be in SI units and can be easily converted to what unit you want.

Parameters:

• x, y1D-arrays

  The x and y coordinates of the grid points

• data1D-array

  The potential field at the grid points

• shapetuple = (nx, ny)

  The shape of the grid

• orderint

  The order of the derivative

• methodstring

  The method used to calculate the derivatives. Options are: 'fft' for the Fast Fourier Transform.

Returns:

• deriv1D-array

  The derivative

geoist.pfm.pftrans.power_density_spectra(x, y, data, shape)

Calculates the Power Density Spectra of a 2D gridded potential field through the FFT:

\[\Phi_{\Delta T}(k_x, k_y) = | F\left{\Delta T \right}(k_x, k_y) |^2\]

注解

Requires gridded data.

注解

x, y, z and height should be in meters.

Parameters:

• x, y1D-arrays

  The x and y coordinates of the grid points

• data1D-array

  The potential field at the grid points

• shapetuple = (nx, ny)

  The shape of the grid

Returns:

• kx, ky2D-arrays

  The wavenumbers of each Power Density Spectra point

• pds2D-array

  The Power Density Spectra of the data

geoist.pfm.pftrans.radial_average_spectrum(kx, ky, pds, max_radius=None, ring_width=None)

Calculates the average of the Power Density Spectra points that falls inside concentric rings built around the origin of the wavenumber coordinate system with constant width.

The width of the rings and the inner radius of the biggest ring can be changed by setting the optional parameters ring_width and max_radius, respectively.

注解

To calculate the radially averaged power density spectra use the outputs of the function power_density_spectra as input of this one.

Parameters:

• kx, ky2D-arrays

  The wavenumbers arrays in the x and y directions

• pds2D-array

  The Power Density Spectra

• max_radiusfloat (optional)

  Inner radius of the biggest ring. By default it's set as the minimum of kx.max() and ky.max(). Making it smaller leaves points outside of the averaging, and making it bigger includes points nearer to the boundaries.

• ring_widthfloat (optional)

  Width of the rings. By default it's set as the largest value of \(\Delta k_x\) and \(\Delta k_y\), being them the equidistances of the kx and ky arrays. Making it bigger gives more populated averages, and making it smaller lowers the ammount of points per ring (use it carefully).

Returns:

• k_radial1D-array

  Wavenumbers of each Radially Averaged Power Spectrum point. Also, the inner radius of the rings.

• pds_radial1D array

  Radially Averaged Power Spectrum

geoist.pfm.pftrans.reduce_to_pole(x, y, data, shape, inc, dec, sinc, sdec)

Reduce total field magnetic anomaly data to the pole.

The reduction to the pole if a phase transformation that can be applied to total field magnetic anomaly data. It "simulates" how the data would be if both the Geomagnetic field and the magnetization of the source were vertical (\(90^\circ\) inclination) (Blakely, 1996).

This functions performs the reduction in the frequency domain (using the FFT). The transform filter is (in the frequency domain):

\[RTP(k_x, k_y) = \frac{|k|}{ a_1 k_x^2 + a_2 k_y^2 + a_3 k_x k_y + i|k|(b_1 k_x + b_2 k_y)}\]

in which \(k_x\) and \(k_y\) are the wave-numbers in the x and y directions and

\[\begin{split}|k| = \sqrt{k_x^2 + k_y^2} \\ a_1 = m_z f_z - m_x f_x \\ a_2 = m_z f_z - m_y f_y \\ a_3 = -m_y f_x - m_x f_y \\ b_1 = m_x f_z + m_z f_x \\ b_2 = m_y f_z + m_z f_y\end{split}\]

\(\mathbf{m} = (m_x, m_y, m_z)\) is the unit-vector of the total magnetization of the source and \(\mathbf{f} = (f_x, f_y, f_z)\) is the unit-vector of the Geomagnetic field.

注解

Requires gridded data.

警告

The magnetization direction of the anomaly source is crucial to the reduction-to-the-pole. Wrong values of *sinc* and *sdec* will lead to a wrong reduction.

Parameters:

• x, y1d-arrays

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

• data1d-array

  The total field anomaly data at each point.

• shapetuple = (nx, ny)

  The shape of the data grid

• inc, decfloats

  The inclination and declination of the inducing Geomagnetic field

• sinc, sdecfloats

  The inclination and declination of the total magnetization of the anomaly source. The total magnetization is the vector sum of the induced and remanent magnetization. If there is only induced magnetization, use the inc and dec of the Geomagnetic field.

Returns:

• rtp1d-array

  The data reduced to the pole.

References:

Blakely, R. J. (1996), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.pftrans.tga(x, y, data, shape, method='fd')

Calculate the total gradient amplitude (TGA).

This the same as the 3D analytic signal of Roest et al. (1992), but we prefer the newer, more descriptive nomenclature suggested by Reid (2012).

The TGA is defined as the amplitude of the gradient vector of a potential field \(T\) (e.g. the magnetic total field anomaly):

\[TGA = \sqrt{ \left(\frac{\partial T}{\partial x}\right)^2 + \left(\frac{\partial T}{\partial y}\right)^2 + \left(\frac{\partial T}{\partial z}\right)^2 }\]

注解

Requires gridded data.

警告

If the data is not in SI units, the derivatives will be in strange units and so will the total gradient amplitude! I strongly recommend converting the data to SI before calculating the TGA is you need the gradient in Eotvos (use one of the unit conversion functions of geoist.utils).

Parameters:

• x, y1D-arrays

  The x and y coordinates of the grid points

• data1D-array

  The potential field at the grid points

• shapetuple = (nx, ny)

  The shape of the grid

• methodstring

  The method used to calculate the horizontal derivatives. Options are: 'fd' for finite-difference (more stable) or 'fft' for the Fast Fourier Transform. The z derivative is always calculated by FFT.

Returns:

• tga1D-array

  The amplitude of the total gradient

References:

Reid, A. (2012), Forgotten truths, myths and sacred cows of Potential Fields Geophysics - II, in SEG Technical Program Expanded Abstracts 2012, pp. 1-3, Society of Exploration Geophysicists.

Roest, W., J. Verhoef, and M. Pilkington (1992), Magnetic interpretation using the 3-D analytic signal, GEOPHYSICS, 57(1), 116-125, doi:10.1190/1.1443174.

geoist.pfm.pftrans.tilt(x, y, data, shape, xderiv=None, yderiv=None, zderiv=None)

Calculates the potential field tilt, as defined by Miller and Singh (1994)

\[tilt(f) = tan^{-1}\left( \frac{ \frac{\partial T}{\partial z}}{ \sqrt{\frac{\partial T}{\partial x}^2 + \frac{\partial T}{\partial y}^2}} \right)\]

When used on magnetic total field anomaly data, works best if the data is reduced to the pole.

It's useful to plot the zero contour line of the tilt to represent possible outlines of the source bodies. Use matplotlib's pyplot.contour or pyplot.tricontour for this.

注解

Requires gridded data if xderiv, yderiv and zderiv are not given.

Parameters:

• x, y1D-arrays

  The x and y coordinates of the grid points

• data1D-array

  The potential field at the grid points

• shapetuple = (nx, ny)

  The shape of the grid. Ignored if xderiv, yderiv and zderiv are given.

• xderiv1D-array or None

  Optional. Values of the derivative in the x direction. If None, will calculated using the default options of derivx

• yderiv1D-array or None

  Optional. Values of the derivative in the y direction. If None, will calculated using the default options of derivy

• zderiv1D-array or None

  Optional. Values of the derivative in the z direction. If None, will calculated using the default options of derivz

Returns:

• tilt1D-array

  The tilt angle of the total field in radians.

References:

Miller, Hugh G, and Vijay Singh. 1994. "Potential Field Tilt --- a New Concept for Location of Potential Field Sources." Journal of Applied Geophysics 32 (2--3): 213-17. doi:10.1016/0926-9851(94)90022-1.

geoist.pfm.pftrans.upcontinue(x, y, data, shape, height)

Upward continuation of potential field data.

Calculates the continuation through the Fast Fourier Transform in the wavenumber domain (Blakely, 1996):

\[F\{h_{up}\} = F\{h\} e^{-\Delta z |k|}\]

and then transformed back to the space domain. \(h_{up}\) is the upward continue data, \(\Delta z\) is the height increase, \(F\) denotes the Fourier Transform, and \(|k|\) is the wavenumber modulus.

注解

Requires gridded data.

注解

x, y, z and height should be in meters.

注解

It is not possible to get the FFT of a masked grid. The default geoist.gridder.interp call using minimum curvature will not be suitable. Use extrapolate=True or algorithm='nearest' to get an unmasked grid.

Parameters:

• x, y1D-arrays

  The x and y coordinates of the grid points

• data1D-array

  The potential field at the grid points

• shapetuple = (nx, ny)

  The shape of the grid

• heightfloat

  The height increase (delta z) in meters.

Returns:

• contarray

  The upward continued data

References:

Blakely, R. J. (1996), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.polyprism module

The potential fields of a homogeneous 3D prism with polygonal cross-section.

geoist.pfm.polyprism.bx(xp, yp, zp, prisms)

x component of magnetic induction of a polygonal prism.

注解

Input units are SI. Output is in nT

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the total field anomaly. Prisms without the physical property 'magnetization' will be ignored. The 'magnetization' must be a vector.

Returns:

• bx: array

  The x component of the magnetic induction

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.by(xp, yp, zp, prisms)

y component of magnetic induction of a polygonal prism.

注解

Input units are SI. Output is in nT

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the total field anomaly. Prisms without the physical property 'magnetization' will be ignored. The 'magnetization' must be a vector.

Returns:

• by: array

  The y component of the magnetic induction

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.bz(xp, yp, zp, prisms)

z component of magnetic induction of a polygonal prism.

注解

Input units are SI. Output is in nT

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the total field anomaly. Prisms without the physical property 'magnetization' will be ignored. The 'magnetization' must be a vector.

Returns:

• bz: array

  The z component of the magnetic induction

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.gxx(xp, yp, zp, prisms)

xx component of the gravity gradient tensor of a polygonal prism.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input values in SI units and output in Eotvos!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the field. Prisms must have the physical property 'density' will be ignored.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.gxy(xp, yp, zp, prisms)

xy component of the gravity gradient tensor of a polygonal prism.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input values in SI units and output in Eotvos!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the field. Prisms must have the physical property 'density' will be ignored.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.gxz(xp, yp, zp, prisms)

xz component of the gravity gradient tensor of a polygonal prism.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input values in SI units and output in Eotvos!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the field. Prisms must have the physical property 'density' will be ignored.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.gyy(xp, yp, zp, prisms)

yy component of the gravity gradient tensor of a polygonal prism.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input values in SI units and output in Eotvos!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the field. Prisms must have the physical property 'density' will be ignored.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.gyz(xp, yp, zp, prisms)

yz component of the gravity gradient tensor of a polygonal prism.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input values in SI units and output in Eotvos!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the field. Prisms must have the physical property 'density' will be ignored.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.gz(xp, yp, zp, prisms)

z component of gravitational acceleration of a polygonal prism.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input values in SI units and output in mGal!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the field. Prisms must have the physical property 'density' will be ignored.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.gzz(xp, yp, zp, prisms)

zz component of the gravity gradient tensor of a polygonal prism.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input values in SI units and output in Eotvos!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismslist of geoist.inversion.mesher.PolygonalPrism

  The model used to calculate the field. Prisms must have the physical property 'density' will be ignored.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.kernelxx(xp, yp, zp, prism)

The xx second-derivative of the kernel function \(\phi\).

\[\phi(x,y,z) = \iiint_\Omega \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

in which

\[r = \sqrt{(x - \nu)^2 + (y - \eta)^2 + (z - \zeta)^2}.\]

This function is used to calculate the gravity gradient tensor, magnetic induction, and total field magnetic anomaly.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input and output values in SI!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.inversion.mesher.PolygonalPrism

  The model used as the integration domain \(\Omega\) of the kernel function.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.kernelxy(xp, yp, zp, prism)

The xy second-derivative of the kernel function \(\phi\).

\[\phi(x,y,z) = \iiint_\Omega \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

in which

\[r = \sqrt{(x - \nu)^2 + (y - \eta)^2 + (z - \zeta)^2}.\]

This function is used to calculate the gravity gradient tensor, magnetic induction, and total field magnetic anomaly.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input and output values in SI!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.inversion.mesher.PolygonalPrism

  The model used as the integration domain \(\Omega\) of the kernel function.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.kernelxz(xp, yp, zp, prism)

The xz second-derivative of the kernel function \(\phi\).

\[\phi(x,y,z) = \iiint_\Omega \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

in which

\[r = \sqrt{(x - \nu)^2 + (y - \eta)^2 + (z - \zeta)^2}.\]

This function is used to calculate the gravity gradient tensor, magnetic induction, and total field magnetic anomaly.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input and output values in SI!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.inversion.mesher.PolygonalPrism

  The model used as the integration domain \(\Omega\) of the kernel function.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.kernelyy(xp, yp, zp, prism)

The yy second-derivative of the kernel function \(\phi\).

\[\phi(x,y,z) = \iiint_\Omega \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

in which

\[r = \sqrt{(x - \nu)^2 + (y - \eta)^2 + (z - \zeta)^2}.\]

This function is used to calculate the gravity gradient tensor, magnetic induction, and total field magnetic anomaly.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input and output values in SI!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.inversion.mesher.PolygonalPrism

  The model used as the integration domain \(\Omega\) of the kernel function.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.kernelyz(xp, yp, zp, prism)

The yz second-derivative of the kernel function \(\phi\).

\[\phi(x,y,z) = \iiint_\Omega \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

in which

\[r = \sqrt{(x - \nu)^2 + (y - \eta)^2 + (z - \zeta)^2}.\]

This function is used to calculate the gravity gradient tensor, magnetic induction, and total field magnetic anomaly.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input and output values in SI!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.inversion.mesher.PolygonalPrism

  The model used as the integration domain \(\Omega\) of the kernel function.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.kernelzz(xp, yp, zp, prism)

The zz second-derivative of the kernel function \(\phi\).

\[\phi(x,y,z) = \iiint_\Omega \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

in which

\[r = \sqrt{(x - \nu)^2 + (y - \eta)^2 + (z - \zeta)^2}.\]

This function is used to calculate the gravity gradient tensor, magnetic induction, and total field magnetic anomaly.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

All input and output values in SI!

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.inversion.mesher.PolygonalPrism

  The model used as the integration domain \(\Omega\) of the kernel function.

Returns:

• resarray

  The effect calculated on the computation points.

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.polyprism.tf(xp, yp, zp, prisms, inc, dec, pmag=None)

The total-field magnetic anomaly of polygonal prisms.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

注解

Input units are SI. Output is in nT

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of geoist.inversion.geometry.PolygonalPrism

  The model used to calculate the total field anomaly. Prisms without the physical property 'magnetization' will be ignored.

• incfloat

  The inclination of the regional field (in degrees)

• decfloat

  The declination of the regional field (in degrees)

• pmag[mx, my, mz] or None

  A magnetization vector. If not None, will use this value instead of the 'magnetization' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

References:

Plouff, D. , 1976, Gravity and magnetic fields of polygonal prisms and applications to magnetic terrain corrections, Geophysics, 41(4), 727-741, doi:10.1190/1.1440645.

geoist.pfm.prism module

Calculate the potential fields of the 3D right rectangular prism.

注解

All input units are SI. Output is in conventional units: SI for the gravitatonal potential, mGal for gravity, Eotvos for gravity gradients, nT for magnetic total field anomalies.

注解

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

Gravity

The gravitational fields are calculated using the formula of Nagy et al. (2000). Available functions are:

• potential

• gx

• gy

• gz

• gxx

• gxy

• gxz

• gyy

• gyz

• gzz

警告

The gxy, gxz, and gyz components have singularities when the computation point is aligned with the corners of the prism on the bottom, east, and north sides, respectively. In these cases, the above functions will move the computation point slightly to avoid these singularities. Unfortunately, this means that the result will not be as accurate on those points.

Magnetic

Available fields are the total-field anomaly (using the formula of Bhattacharyya, 1964) and x, y, z components of the magnetic induction:

• tf

• bx

• by

• bz

Auxiliary Functions

Calculates the second derivatives of the function

\[\phi(x,y,z) = \int\int\int \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

with respect to the variables \(x\), \(y\), and \(z\). In this equation,

\[r = \sqrt{(x - \nu)^2 + (y - \eta)^2 + (z - \zeta)^2}\]

and \(\nu\), \(\eta\), \(\zeta\) are the Cartesian coordinates of an element inside the volume of a 3D prism. These second derivatives are used to calculate the total field anomaly and the gravity gradient tensor components.

• kernelxx

• kernelxy

• kernelxz

• kernelyy

• kernelyz

• kernelzz

References

Bhattacharyya, B. K. (1964), Magnetic anomalies due to prism-shaped bodies with arbitrary polarization, Geophysics, 29(4), 517, doi: 10.1190/1.1439386.

Nagy, D., G. Papp, and J. Benedek (2000), The gravitational potential and its derivatives for the prism: Journal of Geodesy, 74, 552--560, doi: 10.1007/s001900000116.

---

geoist.pfm.prism.bx(xp, yp, zp, prisms, pmag=None)

Calculates the x component of the magnetic induction produced by rectangular prisms.

注解

Input units are SI. Output is in nT

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• prismslist of geoist.mesher.Prism

  The model used to calculate the total field anomaly. Prisms without the physical property 'magnetization' will be ignored. The 'magnetization' must be a vector.

• pmag[mx, my, mz] or None

  A magnetization vector. If not None, will use this value instead of the 'magnetization' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• bx: array

  The x component of the magnetic induction

geoist.pfm.prism.by(xp, yp, zp, prisms, pmag=None)

Calculates the y component of the magnetic induction produced by rectangular prisms.

注解

Input units are SI. Output is in nT

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• prismslist of geoist.mesher.Prism

  The model used to calculate the total field anomaly. Prisms without the physical property 'magnetization' will be ignored. The 'magnetization' must be a vector.

• pmag[mx, my, mz] or None

  A magnetization vector. If not None, will use this value instead of the 'magnetization' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• by: array

  The y component of the magnetic induction

geoist.pfm.prism.bz(xp, yp, zp, prisms, pmag=None)

Calculates the z component of the magnetic induction produced by rectangular prisms.

注解

Input units are SI. Output is in nT

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• prismslist of geoist.mesher.Prism

  The model used to calculate the total field anomaly. Prisms without the physical property 'magnetization' will be ignored. The 'magnetization' must be a vector.

• pmag[mx, my, mz] or None

  A magnetization vector. If not None, will use this value instead of the 'magnetization' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• bz: array

  The z component of the magnetic induction

geoist.pfm.prism.gx(xp, yp, zp, prisms, dens=None)

Calculates the \(g_x\) gravity acceleration component.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in mGal!

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.gxx(xp, yp, zp, prisms, dens=None)

Calculates the \(g_{xx}\) gravity gradient tensor component.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in Eotvos!

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.gxy(xp, yp, zp, prisms, dens=None)

Calculates the \(g_{xy}\) gravity gradient tensor component.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in Eotvos!

警告

This component has singularities when the computation point is aligned with the corners of the prism on the bottom side. In these cases, the computation point slightly to avoid these singularities. Unfortunately, this means that the result will not be as accurate on those points.

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.gxz(xp, yp, zp, prisms, dens=None)

Calculates the \(g_{xz}\) gravity gradient tensor component.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in Eotvos!

警告

This component has singularities when the computation point is aligned with the corners of the prism on the east side. In these cases, the computation point slightly to avoid these singularities. Unfortunately, this means that the result will not be as accurate on those points.

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.gy(xp, yp, zp, prisms, dens=None)

Calculates the \(g_y\) gravity acceleration component.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in mGal!

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.gyy(xp, yp, zp, prisms, dens=None)

Calculates the \(g_{yy}\) gravity gradient tensor component.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in Eotvos!

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.gyz(xp, yp, zp, prisms, dens=None)

Calculates the \(g_{yz}\) gravity gradient tensor component.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in Eotvos!

警告

This component has singularities when the computation point is aligned with the corners of the prism on the north side. In these cases, the computation point slightly to avoid these singularities. Unfortunately, this means that the result will not be as accurate on those points.

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.gz(xp, yp, zp, prisms, dens=None)

Calculates the \(g_z\) gravity acceleration component.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in mGal!

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.gz_ckernel(xp, yp, zp, prisms, dim=None, dens=None, wavelet='db1', levels=1, thresh=0.0001)

Calculates the compressed kernel matrix of gz.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in mGal!

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• rescsr_matrix

  The compressed kernel calculated on xp, yp, zp of prisms.

geoist.pfm.prism.gz_kernel(xp, yp, zp, prisms, dens=None)

Calculates the kernel of gz.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in mGal!

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The kernel calculated on xp, yp, zp of prisms.

geoist.pfm.prism.gzz(xp, yp, zp, prisms, dens=None)

Calculates the \(g_{zz}\) gravity gradient tensor component.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input values in SI units(!) and output in Eotvos!

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.kernelxx(xp, yp, zp, prism)

Calculates the xx derivative of the function

\[\phi(x,y,z) = \int\int\int \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.mesher.Prism

  The model used to calculate the field.

Returns:

• resarray

  The effect calculated on the computation points.

geoist.pfm.prism.kernelxy(xp, yp, zp, prism)

Calculates the xy derivative of the function

\[\phi(x,y,z) = \int\int\int \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.mesher.Prism

  The model used to calculate the field.

Returns:

• resarray

  The effect calculated on the computation points.

geoist.pfm.prism.kernelxz(xp, yp, zp, prism)

Calculates the xz derivative of the function

\[\phi(x,y,z) = \int\int\int \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.mesher.Prism

  The model used to calculate the field.

Returns:

• resarray

  The effect calculated on the computation points.

geoist.pfm.prism.kernelyy(xp, yp, zp, prism)

Calculates the yy derivative of the function

\[\phi(x,y,z) = \int\int\int \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.mesher.Prism

  The model used to calculate the field.

Returns:

• resarray

  The effect calculated on the computation points.

geoist.pfm.prism.kernelyz(xp, yp, zp, prism)

Calculates the yz derivative of the function

\[\phi(x,y,z) = \int\int\int \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.mesher.Prism

  The model used to calculate the field.

Returns:

• resarray

  The effect calculated on the computation points.

geoist.pfm.prism.kernelzz(xp, yp, zp, prism)

Calculates the zz derivative of the function

\[\phi(x,y,z) = \int\int\int \frac{1}{r} \mathrm{d}\nu \mathrm{d}\eta \mathrm{d}\zeta\]

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates of the computation points.

• prismsobject of geoist.mesher.Prism

  The model used to calculate the field.

Returns:

• resarray

  The effect calculated on the computation points.

geoist.pfm.prism.potential(xp, yp, zp, prisms, dens=None)

Calculates the gravitational potential.

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> DOWN.

注解

All input and output values in SI units(!)!

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The density model used to calculate the gravitational effect. Prisms must have the property 'density'. Prisms that don't have this property will be ignored in the computations. Elements of prisms that are None will also be ignored. prisms can also be a PrismMesh.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.prism.tf(xp, yp, zp, prisms, inc, dec, pmag=None)

Calculate the total-field magnetic anomaly of prisms.

注解

Input units are SI. Output is in nT

注解

The coordinate system of the input parameters is to be x -> North, y -> East and z -> Down.

Parameters:

• xp, yp, zparrays

  Arrays with the x, y, and z coordinates of the computation points.

• prismslist of Prism

  The model used to calculate the total field anomaly. Prisms without the physical property 'magnetization' will be ignored. prisms can also be a PrismMesh.

• incfloat

  The inclination of the regional field (in degrees)

• decfloat

  The declination of the regional field (in degrees)

• pmag[mx, my, mz] or None

  A magnetization vector. If not None, will use this value instead of the 'magnetization' property of the prisms. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

geoist.pfm.sphere module

The potential fields of a homogeneous sphere.

geoist.pfm.sphere.bx(xp, yp, zp, spheres, pmag=None)

The x component of the magnetic induction.

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

Input units should be SI. Output is in nT.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the physical property 'magnetization'. Spheres that are None or without 'magnetization' will be ignored. The magnetization is the total (remanent + induced + any demagnetization effects) magnetization given as a 3-component vector.

• pmag[mx, my, mz] or None

  A magnetization vector. If not None, will use this value instead of the 'magnetization' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• bx: array

  The x component of the magnetic induction

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.by(xp, yp, zp, spheres, pmag=None)

The y component of the magnetic induction.

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

Input units should be SI. Output is in nT.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the physical property 'magnetization'. Spheres that are None or without 'magnetization' will be ignored. The magnetization is the total (remanent + induced + any demagnetization effects) magnetization given as a 3-component vector.

• pmag[mx, my, mz] or None

  A magnetization vector. If not None, will use this value instead of the 'magnetization' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• by: array

  The y component of the magnetic induction

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.bz(xp, yp, zp, spheres, pmag=None)

The z component of the magnetic induction.

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

Input units should be SI. Output is in nT.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the physical property 'magnetization'. Spheres that are None or without 'magnetization' will be ignored. The magnetization is the total (remanent + induced + any demagnetization effects) magnetization given as a 3-component vector.

• pmag[mx, my, mz] or None

  A magnetization vector. If not None, will use this value instead of the 'magnetization' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• bzarray

  The z component of the magnetic induction

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.gxx(xp, yp, zp, spheres, dens=None)

The \(g_{xx}\) gravity gradient component.

\[g_{xx}(x, y, z) = \rho 4 \pi \dfrac{radius^3}{3} \dfrac{3 (x - x')^2 - r^2}{r^5}\]

in which \(\rho\) is the density and \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in Eotvos.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the field will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the property 'density'. The ones that are None or without a density will be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.gxy(xp, yp, zp, spheres, dens=None)

The \(g_{xy}\) gravity gradient component.

\[g_{xy}(x, y, z) = \rho 4 \pi \dfrac{radius^3}{3} \dfrac{3(x - x')(y - y')}{r^5}\]

in which \(\rho\) is the density and \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in Eotvos.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the field will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the property 'density'. The ones that are None or without a density will be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.gxz(xp, yp, zp, spheres, dens=None)

The \(g_{xz}\) gravity gradient component.

\[g_{xz}(x, y, z) = \rho 4 \pi \dfrac{radius^3}{3} \dfrac{3(x - x')(z - z')}{r^5}\]

in which \(\rho\) is the density and \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in Eotvos.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the field will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the property 'density'. The ones that are None or without a density will be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.gyy(xp, yp, zp, spheres, dens=None)

The \(g_{yy}\) gravity gradient component.

\[g_{yy}(x, y, z) = \rho 4 \pi \dfrac{radius^3}{3} \dfrac{3(y - y')^2 - r^2}{r^5}\]

in which \(\rho\) is the density and \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in Eotvos.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the field will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the property 'density'. The ones that are None or without a density will be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.gyz(xp, yp, zp, spheres, dens=None)

The \(g_{yz}\) gravity gradient component.

\[g_{yz}(x, y, z) = \rho 4 \pi \dfrac{radius^3}{3} \dfrac{3(y - y')(z - z')}{r^5}\]

in which \(\rho\) is the density and \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in Eotvos.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the field will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the property 'density'. The ones that are None or without a density will be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.gz(xp, yp, zp, spheres, dens=None)

The \(g_z\) gravitational acceleration component.

\[g_z(x, y, z) = \rho 4 \pi \dfrac{radius^3}{3} \dfrac{z - z'}{r^3}\]

in which \(\rho\) is the density and \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in mGal.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the field will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the property 'density'. The ones that are None or without a density will be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.gzz(xp, yp, zp, spheres, dens=None)

The \(g_{zz}\) gravity gradient component.

\[g_{zz}(x, y, z) = \rho 4 \pi \dfrac{radius^3}{3} \dfrac{3(z - z')^2 - r^2}{r^5}\]

in which \(\rho\) is the density and \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in Eotvos.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the field will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the property 'density'. The ones that are None or without a density will be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• resarray

  The field calculated on xp, yp, zp

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.sphere.kernelxx(xp, yp, zp, sphere)

The second x derivative of the kernel function

\[\phi(x,y,z) = \frac{4}{3} \pi radius^3 \frac{1}{r}\]

where \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in SI.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the function will be calculated

• spherefatiando.mesher.Sphere

  The sphere.

Returns:

• resarray

  The function calculated on xp, yp, zp

geoist.pfm.sphere.kernelxy(xp, yp, zp, sphere)

The xy derivative of the kernel function

\[\phi(x,y,z) = \frac{4}{3} \pi radius^3 \frac{1}{r}\]

where \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in SI.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the function will be calculated

• spherefatiando.mesher.Sphere

  The sphere.

Returns:

• resarray

  The function calculated on xp, yp, zp

geoist.pfm.sphere.kernelxz(xp, yp, zp, sphere)

The xz derivative of the kernel function

\[\phi(x,y,z) = \frac{4}{3} \pi radius^3 \frac{1}{r}\]

where \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in SI.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the function will be calculated

• spherefatiando.mesher.Sphere

  The sphere.

Returns:

• resarray

  The function calculated on xp, yp, zp

geoist.pfm.sphere.kernelyy(xp, yp, zp, sphere)

The second y derivative of the kernel function

\[\phi(x,y,z) = \frac{4}{3} \pi radius^3 \frac{1}{r}\]

where \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in SI.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the function will be calculated

• spherefatiando.mesher.Sphere

  The sphere.

Returns:

• resarray

  The function calculated on xp, yp, zp

geoist.pfm.sphere.kernelyz(xp, yp, zp, sphere)

The yz derivative of the kernel function

\[\phi(x,y,z) = \frac{4}{3} \pi radius^3 \frac{1}{r}\]

where \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in SI.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the function will be calculated

• spherefatiando.mesher.Sphere

  The sphere.

Returns:

• resarray

  The function calculated on xp, yp, zp

geoist.pfm.sphere.kernelzz(xp, yp, zp, sphere)

The second z derivative of the kernel function

\[\phi(x,y,z) = \frac{4}{3} \pi radius^3 \frac{1}{r}\]

where \(r = \sqrt{(x - x')^2 + (y - y')^2 + (z - z')^2}\).

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

All input values should be in SI and output is in SI.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the function will be calculated

• spherefatiando.mesher.Sphere

  The sphere.

Returns:

• resarray

  The function calculated on xp, yp, zp

geoist.pfm.sphere.tf(xp, yp, zp, spheres, inc, dec, pmag=None)

The total-field magnetic anomaly.

The anomaly is defined as (Blakely, 1995):

\[\Delta T = |\mathbf{T}| - |\mathbf{F}|,\]

where \(\mathbf{T}\) is the measured field and \(\mathbf{F}\) is a reference (regional) field.

The anomaly of a homogeneous sphere can be calculated as:

\[\Delta T \approx \hat{\mathbf{F}}\cdot\mathbf{B}.\]

where \(\mathbf{B}\) is the magnetic induction produced by the sphere.

The coordinate system of the input parameters is x -> North, y -> East and z -> Down.

Input units should be SI. Output is in nT.

Parameters:

• xp, yp, zparrays

  The x, y, and z coordinates where the anomaly will be calculated

• sphereslist of fatiando.mesher.Sphere

  The spheres. Spheres must have the physical property 'magnetization'. Spheres that are None or without 'magnetization' will be ignored. The magnetization is the total (remanent + induced + any demagnetization effects) magnetization given as a 3-component vector.

• inc, decfloats

  The inclination and declination of the regional field (in degrees)

• pmag[mx, my, mz] or None

  A magnetization vector. If not None, will use this value instead of the 'magnetization' property of the spheres. Use this, e.g., for sensitivity matrix building.

Returns:

• tfarray

  The total-field anomaly

References:

Blakely, R. J. (1995), Potential Theory in Gravity and Magnetic Applications, Cambridge University Press.

geoist.pfm.talwani module

geoist.pfm.tensor module

Utilities for operating on the gradient tensor of potential fields.

Functions

• invariants: Calculates the first (\(I_1\)), second (\(I_2\)), and dimensionless (\(I\)) invariants

• eigen: Calculates the eigenvalues and eigenvectors of an array of gradient tensor measurements

• center_of_mass: Estimate the center of mass of sources from the first eigenvector using the method of Beiki and Pedersen (2010)

Theory

Following Pedersen and Rasmussen (1990), the characteristic polynomial of the gravity gradient tensor \(\mathbf{\Gamma}\) is

\[\lambda^3 + I_1\lambda - I_2 = 0\]

where \(\lambda\) is an eigenvalue and \(I_1\) and \(I_2\) are the two invariants. The dimensionless invariant \(I\) is

\[I = -\dfrac{(I_2/2)^2}{(I_1/3)^3}\]

The invariant \(I\) indicates the dimensionality of the source. \(I = 0\) for 2 dimensional bodies and \(I = 1\) for a monopole.

References

Beiki, M., and L. B. Pedersen (2010), Eigenvector analysis of gravity gradient tensor to locate geologic bodies, Geophysics, 75(6), I37, doi:10.1190/1.3484098

Pedersen, L. B., and T. M. Rasmussen (1990), The gradient tensor of potential field anomalies: Some implications on data collection and data processing of maps, Geophysics, 55(12), 1558, doi:10.1190/1.1442807

---

geoist.pfm.tensor.center_of_mass(x, y, z, eigvec1, windows=1, wcenter=None, wmin=None, wmax=None)

Estimates the center of mass of a source using the 1st eigenvector

Uses the method of Beiki and Pedersen (2010) with an expanding window scheme to get the best estimate and deal with multiple sources.

Parameters:

• x, y, zarrays

  The x, y, and z coordinates of the observation points

• eigvec1array (shape = (N, 3) where N is the number of observations)

  The first eigenvector of the gravity gradient tensor at each observation point

• windowsint

  The number of expanding windows to use

• wcenterlist = [x, y]

  The [x, y] coordinates of the center of the expanding windows. Will default to the middle of the data area if None

• wmin, wmaxfloat

  Minimum and maximum size of the expanding windows. Will default to 10% data area and 100% data area, respectively, if None

Returns:

• [xo, yo, zo]floats

  xo, yo, zo are the coordinates of the estimated center of mass

Examples:

Estimate the center of a sphere using some synthetic data:

>>> from geoist import gridder
>>> from geoist.inversion.mesher import Sphere
>>> from geoist.pfm import sphere, tensor
>>> # Generate synthetic data using a sphere model
>>> # The center of the sphere is (-100, 0, 100)
>>> model = [Sphere(-100, 20, 100, 100, {'density':1000})]
>>> x, y, z = gridder.regular((-500, 500, -500, 500), (20, 20), z=-100)
>>> data = [sphere.gxx(x, y, z, model),
...         sphere.gxy(x, y, z, model),
...         sphere.gxz(x, y, z, model),
...         sphere.gyy(x, y, z, model),
...         sphere.gyz(x, y, z, model),
...         sphere.gzz(x, y, z, model)]
>>> # Get the first eigenvector
>>> eigenvals, eigenvecs = tensor.eigen(data)
>>> # Now estimate the center of mass
>>> cm = tensor.center_of_mass(x, y, z, eigenvecs[0])
>>> print('{:.1f} {:.1f} {:.1f}'.format(cm[0], cm[1], cm[2]))
-100.0 20.0 100.0

geoist.pfm.tensor.eigen(tensor)

Calculates the eigenvalues and eigenvectors of the gradient tensor.

注解

The coordinate system used is x->North, y->East, z->Down

Parameters:

• tensorlist

  A list of arrays with the 6 components of the gradient tensor measured on a set of points. The order of the list should be: [gxx, gxy, gxz, gyy, gyz, gzz]

Returns:

• resultlist of lists

  The eigenvalues and eigenvectors at each observation point. [[eigval1, eigval2, eigval3], [eigvec1, eigvec2, eigvec3]]

  • eigval1,2,3array

    The first, second, and third eigenvalues

  • eigvec1,2,3array (shape = (N, 3) where N is the number of points)

    The first, second, and third eigenvectors

Example:

>>> tensor = [[2], [0], [0], [3], [0], [1]]
>>> eigenvals, eigenvecs = eigen(tensor)
>>> print eigenvals[0], eigenvecs[0]
[ 3.] [[ 0.  1.  0.]]
>>> print eigenvals[1], eigenvecs[1]
[ 2.] [[ 1.  0.  0.]]
>>> print eigenvals[2], eigenvecs[2]
[ 1.] [[ 0.  0.  1.]]

geoist.pfm.tensor.invariants(tensor)

Calculates the first, second, and dimensionless invariants of the gradient tensor.

注解

The coordinate system used is x->North, y->East, z->Down

Parameters:

• tensorlist

  A list of arrays with the 6 components of the gradient tensor measured on a set of points. The order of the list should be: [gxx, gxy, gxz, gyy, gyz, gzz]

Returns:

• invariantslist = [\(I_1\), \(I_2\), \(I\)]

  The invariants calculated for each point

geoist.pfm.tesseroid module

Forward model the gravitational fields of a tesseroid (spherical prism).

Functions in this module calculate the gravitational fields of a tesseroid with respect to the local North-oriented coordinate system of the computation point. See the figure below.

A tesseroid in a geocentric coordinate system (X, Y, Z). Point P is a computation point with associated local North-oriented coordinate system (x, y, z). Image by L. Uieda (doi:10.6084/m9.figshare.1495525).

Coordinate systems

The gravitational attraction and gravity gradient tensor are calculated with respect to the local coordinate system of the computation point. This system has x -> North, y -> East, z -> up (radial direction).

警告

The \(g_z\) component is an exception to this. In order to conform with the regular convention of z-axis pointing toward the center of the Earth, this component only is calculated with z -> Down. This way, gravity anomalies of tesseroids with positive density are positive, not negative.

Gravity

Forward modeling of gravitational fields is performed by functions:

potential, gx, gy, gz, gxx, gxy, gxz, gyy, gyz, gzz

The fields are calculated using Gauss-Legendre Quadrature integration and the adaptive discretization algorithm of Uieda et al. (2016). The accuracy of the integration is controlled by the ratio argument. Larger values cause finer discretization and more accuracy but slower computation times. The defaults values are the ones suggested in the paper and guarantee an accuracy of approximately 0.1%.

警告

The integration error may be larger than this if the computation points are closer than 1 km of the tesseroids. This effect is more significant in the gravity gradient components.

引用

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.gx(lon, lat, height, model, dens=None, ratio=1.6, njobs=1, pool=None)

Calculate the North component of the gravitational attraction.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in mGal

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.gxx(lon, lat, height, model, dens=None, ratio=8, njobs=1, pool=None)

Calculate the xx component of the gravity gradient tensor.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in Eotvos

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.gxy(lon, lat, height, model, dens=None, ratio=8, njobs=1, pool=None)

Calculate the xy component of the gravity gradient tensor.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in Eotvos

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.gxz(lon, lat, height, model, dens=None, ratio=8, njobs=1, pool=None)

Calculate the xz component of the gravity gradient tensor.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in Eotvos

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.gy(lon, lat, height, model, dens=None, ratio=1.6, njobs=1, pool=None)

Calculate the East component of the gravitational attraction.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in mGal

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.gyy(lon, lat, height, model, dens=None, ratio=8, njobs=1, pool=None)

Calculate the yy component of the gravity gradient tensor.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in Eotvos

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.gyz(lon, lat, height, model, dens=None, ratio=8, njobs=1, pool=None)

Calculate the yz component of the gravity gradient tensor.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in Eotvos

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.gz(lon, lat, height, model, dens=None, ratio=1.6, njobs=1, pool=None)

Calculate the radial component of the gravitational attraction.

警告

In order to conform with the regular convention of positive density giving positive gz values, this component only is calculated with z axis -> Down.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in mGal

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.gzz(lon, lat, height, model, dens=None, ratio=8, njobs=1, pool=None)

Calculate the zz component of the gravity gradient tensor.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in Eotvos

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.tesseroid.potential(lon, lat, height, model, dens=None, ratio=1, njobs=1, pool=None)

Calculate the gravitational potential due to a tesseroid model.

警告

Tesseroids with dimensions < 10 cm will be ignored to avoid numerical errors.

Implements the method of Uieda et al. (2016).

Parameters:

• lon, lat, heightarrays

  Arrays with the longitude, latitude and height coordinates of the computation points.

• modellist of Tesseroid

  The density model used to calculate the gravitational effect. Tesseroids must have the property 'density'. Those that don't have this property will be ignored in the computations. Elements that are None will also be ignored.

• densfloat or None

  If not None, will use this value instead of the 'density' property of the tesseroids. Use this, e.g., for sensitivity matrix building.

• ratiofloat

  Will divide each tesseroid until the distance between it and the computation points is < ratio*size of tesseroid. Used to guarantee the accuracy of the numerical integration.

• njobsint

  Split the computation into njobs parts and run it in parallel using multiprocessing. If njobs=1 will run the computation in serial.

• poolNone or multiprocessing.Pool object

  If not None, will use this pool to run the computation in parallel instead of creating a new one. You must still specify njobs as the number of processes in the pool. Use this to avoid spawning processes on each call to this functions, which can have significant overhead.

Returns:

• resarray

  The calculated field in SI units

References:

Uieda, L., V. Barbosa, and C. Braitenberg (2016), Tesseroids: Forward-modeling gravitational fields in spherical coordinates, Geophysics, F41-F48, doi:10.1190/geo2015-0204.1

geoist.pfm.test1 module

Created on Wed May 1 18:44:52 2019

@author: chens

geoist.pfm.tide module

Name : tide.py Created on : 2018/11/03 17:00 Author : Steve Chen <chenshi@cea-igp.ac.cn> Affiliation : Institute of Geophysics, CEA. Version : 0.1.0 Copyright : Copyright (C) 2018-2020 GEOIST Development Team. All Rights Reserved. License : Distributed under the MIT License. See LICENSE.txt for more info. Github : https://igp-gravity.github.io/ Description : code for UTC time, 1.16 don't multiply.

class geoist.pfm.tide.TideModel

基类：object

Class to encapsulate the Longman 1959 tide model.

calculate_julian_century(timestamp)

Calculate the julian century and hour.

Take a datetime object and calculate the decimal Julian century and floating point hour. This is in reference to noon on December 31, 1899 as stated in the Longman paper.

参数

timestamp (datetime) -- Time stamp to convert

返回

Julian century and hour

返回类型

float, float

plot()

Plot the model results.

Make a simple plot of the gravitational tide results from the model run.

run_model()

Run the model for a range of times.

Runs the tidal model beginning at start_time with time steps of increment seconds for days.

solve_longman(lat, lon, alt, time)

Solve the tide model.

Given the location and datetime object, computes the current gravitational tide and associated quantities. Latitude and longitude and in the traditional decimal notation, altitude is in meters, time is a datetime object.

参数

• lat (float) -- latitude (in degrees)

• lon (float) -- longitude (in degrees)

• alt (float) -- altitude (in meters)

• time (datetime) -- time at which to solve the model

返回

lunar, solar, and total gravitational tides

返回类型

float, float, float

write(fname)

Write model results to file.

Write results out of a file for later analysis or reading into another method for analysis/correction of data.

参数

fname (string) -- name of file to save

Module contents