Package Documentation

python interface for global geomagnetic field models (pymagglobal)

pymagglobal serves the purpose of replacing some Fortran scripts, which are used in the geomagnetism community to evaluate global field models. It can be applied to all cubic-spline based geomagnetic field models stored in the same file format as gufm1 or the CALSxk model series.

The package exposes several methods that can be applied to models or splines representing a model. To acces the models, use the Model class or its attribute Model.splines. With the models dictionary, several built-in models can be accessed, see also the list of included models.

class pymagglobal.Model(name, fname=None)[source]

Class for representing a global archeomagnetic field model.

Parameters
  • name (str) – The name of the model. If fname is none, the built-in model will be used.

  • fname (str, optional) – The path leading to a model file in the format of gufm1 and the CALSxk- series.

t_min, t_max[source]

The lower and upper bounds for the model time series. Above and below these the time series will be extrapolated and not robust.

Type

float

l_max[source]

The maximal spherical harmonics degree included in the model.

Type

int

knots[source]

The epochs of the spline knots.

Type

array

coeffs[source]

A large array containing coefficients at the (inner) knots.

Type

array

splines[source]

A BSpline instance representing the model.

Type

scipy.interpolate.BSpline

cov_splines[source]

A BSPline instance giving the covariances. Only available for certain models. Default is None.

Type

scipy.interpolate.BSpline

valid_degrees_orders(degrees, orders)[source]

Check whether the given degrees and orders are valid for the model and calculate the corresponding indices.

Parameters
  • degrees (int or array-like of ints) – The degrees.

  • orders (int or array-like of ints) – The orders.

Returns

The standard order indices corresponding to the given degrees and orders.

Return type

list of ints

valid_epoch(epoch)[source]

Check whether the given epoch is in the range of the model.

Parameters

epoch (float) – The epoch in years.

Returns

The epoch in years.

Return type

float

pymagglobal.built_in_models() dict[source]

Returns a dictionary of models for all model-files in the dat-folder, together with the paths to the model-files.

pymagglobal.coefficients(epoch, splines, cov_splines=None, check=True)[source]

Evaluate splines at an epoch and return the coefficients, together with the appropriate degrees and orders.

Parameters
  • epoch (float or array) – The epoch(s) at which to return the coefficients, in years.

  • splines (scipy.interpolate.BSpline or Model) – An instance of Model or splines specifying the model.

  • cov_splines (scipy.interpolate.BSpline, optional) – A BSpline interpolating the covariance matrices of the coefficients. If a model is passed to splines, this argument will be ignored and set, if the model’s covariances are available. If cov_splines is given, the standard deviations are returned as a fourth return value.

  • check (bool, optional) – If a Model is given, check the input validity.

Returns

  • list – The list of coefficient degrees

  • list – The list of coefficient orders

  • ndarray – Coefficients of the model for the given epoch, corresponding to the lists of degrees and orders.

pymagglobal.dipole_series(times, splines, cov_splines=None, check=True)[source]

Create a dipole-moment time series from a splines object.

Parameters
  • times (array-like) – The times for which to create the dipole series.

  • splines (scipy.interpolate.BSpline or Model) – An instance of Model or splines specifying the model.

  • cov_splines (scipy.interpolate.BSpline, optional) – A BSpline interpolating the covariance matrices of the coefficients. If a model is passed to splines, this argument will be ignored and set, if the model’s covariances are available. If cov_splines is given, the 16- and 84-percentiles are returned along with the series.

  • check (bool, optional) – If a Model is given, check the input validity.

Returns

An array containing the dipole-moment time series. If cov_splines is given, the 16- and 84-percentiles are returned along with the series.

Return type

ndarray

pymagglobal.field(z_at, splines, cov_splines=None, field_type='nez', check=True, inp_gd=True, out_gd=True)[source]

Evaluate coefficient splines at locations z_at and return the field.

Parameters
  • z_at (array-like of shape (4, n) in Fortran order) –

    The locations at which to evaluate the splines. Should be as expected by utils.dsh_basis, i.e.

    • z_at[0] contains co-latitudes in degrees.

    • z_at[1] contains longitudes in degrees.

    • z_at[2] contains radii in km.

    • z_at[3] contains dates in years.

    You can use pymagglobal.utils.get_z_at to generate an array of input points.

  • splines (scipy.interpolate.BSpline or Model) – An instance of Model or splines specifying the model.

  • cov_splines (scipy.interpolate.BSpline, optional) – A BSpline interpolating the covariance matrices of the coefficients. If a model is passed to splines, this argument will be ignored and set, if the model’s covariances are available. If cov_splines is given, the standard deviations are returned as a second return value.

  • field_type ({'dif', 'nez'}) – The type of output. May be either ‘nez’ (default) for north, east and down component or ‘dif’ for declination, inclination and intensity.

  • check (bool, optional) – If a Model is given, check the input validity.

  • inp_gd (bool, optional) – Whether the inputs are given in a geodetic reference frame.

  • out_gd (bool, optional) – Whether the outputs should be given in a geodetic reference frame.

Returns

An array containing the field components at the input locations. If cov_splines is given, the standard deviations are returned as a second return value of identical shape.

Return type

array of shape (3, n)

pymagglobal.file2splines(fname)[source]

Read a coefficient file in the gufm1/CAlSxk format and return a scipy.interpolate.BSpline object representing the model.

pymagglobal.local_curve(times, loc, splines, cov_splines=None, field_type='dif', check=True)[source]

Create local curves from a splines object.

Parameters
  • times (array-like) – The times for which to create the local curve.

  • loc (tuple) – lat, lon tuple of the location at which to create the local curve.

  • splines (scipy.interpolate.BSpline or Model) – An instance of Model or splines specifying the model.

  • field_type ({'dif', 'nez'}) – The type of the local curves. May be either ‘dif’ (default) for declination, inclination and intensity or ‘nez’ for north, east, down.

  • cov_splines (scipy.interpolate.BSpline, optional) – A BSpline interpolating the covariance matrices of the coefficients. If a model is passed to splines, this argument will be ignored and set, if the model’s covariances are available.

  • check (bool, optional) – If a Model is given, check the input validity.

Returns

  • ndarray – The first component local curve. Either declination or north, depending on the field_type kwarg. If cov_splines is given, a tuple containing the component and the standard deviation is returned.

  • ndarray – The second component local curve. Either inclination or east, depending on the field_type kwarg. If cov_splines is given, a tuple containing the component and the standard deviation is returned.

  • ndarray – The third component local curve. Either intensity or down, depending on the field_type kwarg. If cov_splines is given, a tuple containing the component and the standard deviation is returned.

pymagglobal.power_spectrum(epoch, splines, check=True)[source]

Evaluate splines at an epoch and calculate the spatial power spectrum.

Parameters
  • epoch (float or array) – The epoch(s) at which to return the coefficients, in years.

  • splines (scipy.interpolate.BSpline or Model) – An instance of Model or splines specifying the model.

  • check (bool, optional) – If a Model is given, check the input validity.

Returns

  • list – The list of coefficient degrees

  • array – The power spectrum over the given coefficients, with one power per degree.

pymagglobal.secular_variation(epoch, splines, cov_splines=None, check=True)[source]

Evaluate splines at an epoch and return the secular variation coefficients, together with the appropriate degrees and orders.

Parameters
  • epoch (float or array) – The epoch(s) at which to return the coefficients, in years.

  • splines (scipy.interpolate.BSpline or Model) – An instance of Model or splines specifying the model.

  • check (bool, optional) – If a Model is given, check the input validity.

Returns

  • list – The list of coefficient degrees

  • list – The list of coefficient orders

  • ndarray – Secular variation coefficients of the model for the given epoch, corresponding to the lists of degrees and orders.

Utilities

Utility functions for pymagglobal.

class pymagglobal.utils.WGS84[source]

A simple object to store WGS84 parameters.

pymagglobal.utils.dsh_basis(lmax, z, out, R=6371.2)[source]

Write the magnetic field basis functions, evaluated at points z, into the array out. These are basically the derivatives of the spherical harmonics in spherical coordinates with some scaling factors applied. This implementation is based on a specific recursion of the Legendre polynomials, to guarantee sane behavior at the poles. See [Du] for further details.

This function does not return anything. Its interface is chosen to be replacable by a C-based function from the related pyfield library at a later stage [FieldTools].

Parameters
  • lmax (int) – The maximal spherical harmonics degree.

  • z (array of shape (3, n)) –

    The points at which to evaluate the basis functions, given as
    • z[0, :] contains the co-latitude in degrees.

    • z[1, :] contains the longitude in degrees.

    • z[2, :] contains the radius in kilometers.

  • out (array of shape (lmax2N(lmax), 3*n)) –

    The output array in which the basis functions are stored, as follows:
    • out[i, 0::3] contains the North component of the basis

      corresponding to degree l and order m, where i = lm2i(l, m).

    • out[i, 1::3] contains the East component of the basis.

    • out[i, 2::3] contains the Down component of the basis.

References

Du

J. Du, “Non-singular spherical harmonic expressions of geomagnetic vector and gradient tensor fields in the local north-oriented reference frame.”, Geosci. Model Dev., vol. 8, pages 1979-1990, 2014.

FieldTools

H. Matuschek and S. Mauerberger, “FieldTools - A toolbox for manipulating vector fields on the sphere”, GFZ Data Services, 2019. DOI: 10.5880/fidgeo.2019.033

pymagglobal.utils.equi_sph(n, twopi=True)[source]
Regularly place points on the surface of the unit sphere [Deserno]. the

number of output points may differ from n.

Parameters
  • n (int) – The approximate number of points

  • twopi (bool, optional) – Whether to include 2*pi as a duplicate of 0 for the longitutinal angle (useful for plotting purposes).

Returns

polar and azimutal angles in radians of m points, that are equally spaced on the sphere. m is approximately n.

Return type

array of shape (2, m)

References

Deserno

M. Deserno, “How to generate equidistributed points on the surface of a sphere”, Max-Planck-Institut für Polymerforschung, 2004.

pymagglobal.utils.geodetic2geocentric(z, ellipsoid=<pymagglobal.utils.WGS84 object>)[source]

Assuming locations z are given in geodetic coordinates, transform them to geocentric ones.

Details are scattered in this article: https://en.wikipedia.org/wiki/Geographic_coordinate_conversion

Parameters
  • z (array of shape (3, n)) –

    The points at which to evaluate the basis functions, given as
    • z[0, :] contains the geodetic co-latitude in degrees.

    • z[1, :] contains the geodetic longitude in degrees.

    • z[2, :] contains the geodetic altitude in kilometers, wrt. the

      Earth’s center.

  • ellipsoid (object, optional) – An object representing the reference ellipsoid. Should provide ellipsoid parameters a and f. Default is WGS84.

Returns

z_gc – The transformed locations, in geocentric coordinates.

Return type

array of shape (3, n)

pymagglobal.utils.get_z_at(n, R=6371.2, t=1900.0, random=False, twopi=True)[source]
Get input points on the sphere. This is a convenience routine to allow

easy construction of field maps and synthetic data from the models.

Parameters
  • n (int) – The approximate number of points. Due to the points being equally spaced on the sphere, the actual number may be slightly higher.

  • R (float, optional) – The radius of the sphere. By default this is the Earth’s radius.

  • t (float, optional) – The epoch at which the points are generated.

  • random (bool, optional) – If true, exactly n random points are returned. This is useful for generating synthetic data.

  • twopi (bool, optional) – Whether to include 2*pi as a duplicate of 0 for the longitutinal angle (useful for plotting purposes). Only used if random is False.

Returns

z_at

  • z_at[0] contains co-latitudes in degrees.

  • z_at[1] contains longitudes in degrees.

  • z_at[2] contains radii in km.

  • z_at[3] contains dates in years.

n’ is approximately n.

Return type

array of shape (4, n’)

pymagglobal.utils.grad_d(n, e, z)[source]

Calculate the gradient of D.

pymagglobal.utils.grad_i(n, e, z)[source]

Calculate the gradient of I.

pymagglobal.utils.i2lm_l(i)[source]

Returns the degree l of a Gauss coefficient at index i, in the “standard order”. The “standard order” is as follows

i

l

m

0

1

0

1

1

1

2

1

-1

3

2

0

4

2

1

5

2

-1

6

2

2

7

2

-2

Parameters

i (int) – The index of a Gauss coefficient.

Returns

The corresponding degree.

Return type

int

Examples

>>> i2lm_l(0)
1
>>> i2lm_l(21)
4
pymagglobal.utils.i2lm_m(i)[source]

Returns the order m of a Gauss coefficient at index i, in the “standard order”. The “standard order” is as follows

i

l

m

0

1

0

1

1

1

2

1

-1

3

2

0

4

2

1

5

2

-1

6

2

2

7

2

-2

Parameters

i (int) – The index of a Gauss coefficient.

Returns

The corresponding order.

Return type

int

Examples

>>> i2lm_m(0)
0
>>> i2lm_m(21)
-3
pymagglobal.utils.lm2i(ell, m)[source]

Returns the “standard order” index corresponding to a Gauss coefficient of degree l and order m. The “standard order” is as follows

i

l

m

0

1

0

1

1

1

2

1

-1

3

2

0

4

2

1

5

2

-1

6

2

2

7

2

-2

Parameters
  • ell (int) – The degree of a Gauss coefficient.

  • m (int) – The order of a Gauss coefficient.

Returns

The corresponding index.

Return type

int

Examples

>>> lm2i(1, 0)
0
>>> lm2i(4, -1)
17
pymagglobal.utils.lmax2N(lmax)[source]

Returns the number of Gauss coefficients up to and including the degree lmax.

Parameters

lmax (int) – The maximal spherical harmonic degree.

Returns

The number of Gauss coefficients.

Return type

int

pymagglobal.utils.nez2dif(n, e, z)[source]

Transform the magnetic field components north, east and vertically downward to declination, inclination and intensity.

Parameters
  • n – Pointing towards north.

  • e – Pointing towards east.

  • z – Pointing radially downwards.

Returns

  • ndarray – Declination.

  • ndarray – Inclination.

  • ndarray – Intensity.

Example

>>> nez2dif(1, 0, 0)
(0.0, 0.0, 1.0)