chaosmagpy.chaos.BaseModel

class chaosmagpy.chaos.BaseModel(name, breaks=None, order=None, coeffs=None, source=None, meta=None)[source]

Bases: chaosmagpy.chaos.Base

Class for piecewise polynomial spherical harmonic models.

Parameters
namestr

User specified name of the model.

breaksndarray, shape (m+1,)

Break points (m pieces plus one) for the piecewise polynomial representation of the field model in modified Julian date format. If a single break point is given, it is appended to itself to have two points defining the model interval (single piece).

orderint, positive

Order k of the polynomial pieces (e.g. 1 = constant, 4 = cubic).

coeffsndarray, shape (k, m, nmax * (nmax + 2))

Coefficients of the piecewise polynomial representation of the field model.

source{‘internal’, ‘external’}

Internal or external source (defaults to 'internal')

metadict, optional

Dictionary containing additional information about the model if available.

Attributes
breaksndarray, shape (m+1,)

Break points (m pieces plus one) for the piecewise polynomial representation of the field model in modified Julian date format.

piecesint, positive

Number m of intervals given by break points in breaks.

orderint, positive

Order k of the polynomial pieces (e.g. 1 = constant, 4 = cubic).

nmaxint, positive

Maximum spherical harmonic degree of the field model.

dimint, nmax * (nmax + 2)

Dimension of the model.

coeffsndarray, shape (k, m, nmax * (nmax + 2))

Coefficients of the time-dependent field.

source{‘internal’, ‘external’}

Internal or external source (defaults to 'internal')

metadict, optional

Dictionary containing additional information about the model if available.

Methods

from_bspline(name, knots, coeffs, order[, ...])

Return BaseModel instance from a B-spline representation.

from_shc(filepath, *[, name, leap_year, ...])

Return BaseModel instance by loading a model from an shc-file.

plot_maps(time, radius, **kwargs)

Plot global maps of the field components.

plot_power_spectrum(time, **kwargs)

Plot the spatial power spectrum.

plot_timeseries(radius, theta, phi, **kwargs)

Plot the time series of the time-dependent field components at a specific location.

power_spectrum(time[, radius])

Compute the spatial power spectrum.

synth_coeffs(time, *[, nmax, deriv, extrapolate])

Compute the coefficients from the piecewise polynomial representation.

synth_values(time, radius, theta, phi, *[, ...])

Compute magnetic components from the field model.

classmethod from_bspline(name, knots, coeffs, order, source=None, meta=None)[source]

Return BaseModel instance from a B-spline representation.

Parameters
namestr

User specified name of the model.

knotsndarray, shape (N,)

B-spline knots. Knots must have endpoint multiplicity equal to order. Zero-pad coeffs if needed.

coeffsndarray, shape (M, D)

Bspline coefficients for the M B-splines parameterizing D dimensions.

orderint

Order of the B-spline.

source{‘internal’, ‘external’}

Internal or external source (defaults to 'internal')

metadict, optional

Dictionary containing additional information about the model.

Returns
modelBaseModel

Class BaseModel instance.

classmethod from_shc(filepath, *, name=None, leap_year=None, source=None, meta=None)[source]

Return BaseModel instance by loading a model from an shc-file.

Parameters
filepathstr

Path to shc-file.

namestr, optional

User defined name of the model. Defaults to the filename without the file extension.

leap_year{False, True}, optional

Take leap year in time conversion into account. By default set to False, so that a conversion factor of 365.25 days per year is used.

source{‘internal’, ‘external’}

Internal or external source (defaults to 'internal')

metadict, optional

Dictionary containing additional information about the model.

Returns
modelBaseModel

Class BaseModel instance.

plot_maps(time, radius, **kwargs)[source]

Plot global maps of the field components.

Parameters
timendarray, shape (), (1,) or float

Time given as MJD2000 (modified Julian date).

radiusndarray, shape (), (1,) or float

Array containing the radius in kilometers.

Returns
B_radius, B_theta, B_phi

Global map of the radial, colatitude and azimuthal field components.

Other Parameters
nmaxint, positive, optional

Maximum degree of the harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).

derivint, positive, optional

Derivative in time (default is 0). For secular variation, choose deriv=1.

**kwargskeywords

Other options are passed to plot_utils.plot_maps() function.

plot_power_spectrum(time, **kwargs)[source]

Plot the spatial power spectrum.

Parameters
timefloat

Time in modified Julian date.

plot_timeseries(radius, theta, phi, **kwargs)[source]

Plot the time series of the time-dependent field components at a specific location.

Parameters
radiusndarray, shape (), (1,) or float

Radius of station in kilometers.

thetandarray, shape (), (1,) or float

Colatitude in degrees \([0^\circ, 180^\circ]\).

phindarray, shape (), (1,) or float

Longitude in degrees.

Returns
B_radius, B_theta, B_phi

Time series plot of the radial, colatitude and azimuthal field components.

Other Parameters
nmaxint, positive, optional

Maximum degree of the harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).

derivint, positive, optional

Derivative in time (default is 0). For secular variation, choose deriv=1.

extrapolate{‘linear’, ‘spline’, ‘constant’, ‘off’}, optional

Extrapolate to times outside of the model bounds. Defaults to 'linear'.

**kwargskeywords

Other options to pass to plot_utils.plot_timeseries() function.

power_spectrum(time, radius=None, **kwargs)[source]

Compute the spatial power spectrum.

Parameters
timendarray, shape (…)

Time in modified Julian date.

radiusfloat, optional

Radius in kilometers (defaults to mean Earth’s surface defined in basicConfig['r_surf']).

Returns
R_nndarray, shape (…, nmax)

Spatial power spectrum of the spherical harmonic expansion up to degree nmax.

Other Parameters
nmaxint, positive, optional

Maximum degree of the harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).

derivint, positive, optional

Derivative in time (default is 0). For secular variation, choose deriv=1.

**kwargskeywords

Other options to pass to BaseModel.synth_coeffs() method.

synth_coeffs(time, *, nmax=None, deriv=None, extrapolate=None)[source]

Compute the coefficients from the piecewise polynomial representation.

Parameters
timendarray, shape (…) or float

Array containing the time in modified Julian dates.

nmaxint, positive, optional

Maximum degree of the harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).

derivint, positive, optional

Derivative in time (defaults to 0). For secular variation, choose deriv=1.

extrapolate{‘linear’, ‘quadratic’, ‘cubic’, ‘spline’, ‘constant’, ‘off’} or int, optional

Extrapolate to times outside of the model bounds. Specify polynomial degree as string or any order as integer. Defaults to 'linear' (equiv. to order 2 polynomials).

Value

Description

‘constant’

Use degree zero polynomial only (extrapolate=1).

‘linear’

Use degree-1 polynomials (extrapolate=2).

‘quadratic’

Use degree-2 polynomials (extrapolate=3).

‘cubic’

Use degree-3 polynomials (extrapolate=4).

‘spline’

Use all degree polynomials.

‘off’

Return NaN outside model bounds (extrapolate=0).

Returns
coeffsndarray, shape (…, nmax * (nmax + 2))

Array containing the coefficients.

synth_values(time, radius, theta, phi, *, nmax=None, deriv=None, grid=None, extrapolate=None)[source]

Compute magnetic components from the field model.

Parameters
timendarray, shape (…) or float

Array containing the time in modified Julian date.

radiusndarray, shape (…) or float

Radius in kilometers.

thetandarray, shape (…) or float

Colatitude in degrees.

phindarray, shape (…) or float

Longitude in degrees.

nmaxint, positive, optional

Maximum degree of the harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).

derivint, positive, optional

Derivative in time (defaults to 0). For secular variation, choose deriv=1.

gridbool, optional

If True, field components are computed on a regular grid. Arrays theta and phi must have one dimension less than the output grid since the grid will be created as their outer product.

extrapolate{‘linear’, ‘quadratic’, ‘cubic’, ‘spline’, ‘constant’, ‘off’} or int, optional

Extrapolate to times outside of the model bounds. Specify polynomial degree as string or any order as integer. Defaults to 'linear' (equiv. to order 2 polynomials).

Returns
B_radius, B_theta, B_phindarray, shape (…)

Radial, colatitude and azimuthal field components.