chaosmagpy.chaos.CHAOS¶
- class chaosmagpy.chaos.CHAOS(breaks, order=None, *, coeffs_tdep=None, coeffs_static=None, coeffs_sm=None, coeffs_gsm=None, breaks_delta=None, coeffs_delta=None, breaks_euler=None, coeffs_euler=None, breaks_cal=None, coeffs_cal=None, name=None, meta=None)[source]¶
Bases:
object
Class for the time-dependent geomagnetic field model CHAOS.
- Parameters:
- breaksndarray, shape (m+1,)
Break points for piecewise polynomial representation of the time-dependent internal (i.e. large-scale core) field in modified Julian date format.
- orderint, positive
Order k of polynomial pieces (e.g. 4 = cubic) of the time-dependent internal field.
- coeffs_tdepndarray, shape (k, m,
n_tdep
* (n_tdep
+ 2)) Coefficients of the time-dependent internal field as piecewise polynomial.
- coeffs_staticndarray, shape (
n_static
* (n_static
+ 2),) Coefficients of the static internal (i.e. small-scale crustal) field.
- coeffs_smndarray, shape (
n_sm
* (n_sm
+ 2),) Coefficients of the static external field in SM coordinates.
- coeffs_gsmndarray, shape (
n_gsm
* (n_gsm
+ 2),) Coefficients of the static external field in GSM coordinates.
- breaks_deltadict with ndarrays, shape (\(m_q\) + 1,)
Breaks of baseline corrections of static external field in SM coordinates. The dictionary keys are
'q10'
,'q11'
,'s11'
.- coeffs_deltadict with ndarrays, shape (1, \(m_q\))
Coefficients of baseline corrections of static external field in SM coordinates. The dictionary keys are
'q10'
,'q11'
,'s11'
.- breaks_eulerdict with ndarrays, shape (\(m_e\) + 1,)
Dictionary containing satellite name as key and corresponding break vectors of Euler angles (keys are
'oersted'
,'champ'
,'sac_c'
,'swarm_a'
,'swarm_b'
,'swarm_c'
,'cryosat-2_1'
).- coeffs_eulerdict with ndarrays, shape (1, \(m_e\), 3)
Dictionary containing satellite name as key and arrays of the Euler angles alpha, beta and gamma as trailing dimension (keys are
'oersted'
,'champ'
,'sac_c'
,'swarm_a'
,'swarm_b'
,'swarm_c'
,'cryosat-2_1'
).- breaks_caldict with ndarrays, shape (\(m_c\) + 1,)
Dictionary containing satellite name as key and corresponding break vectors for the calibration parameters (keys are
'cryosat-2_1'
).- coeffs_caldict with ndarrays, shape (1, \(m_c\), 3)
Dictionary containing satellite name as key and arrays of the 9 basic calibration parameters (3 offsets, 3 sensitivities, 3 non-orthogonality angles) (keys are
'cryosat-2_1'
).- namestr, optional
User defined name of the model. Defaults to
'CHAOS-<version>'
, where <version> is the version inbasicConfig['params.version']
.- metadict, optional
Dictionary containing additional information about the model.
Examples
Load for example the mat-file
CHAOS-6-x7.mat
in the current working directory like this:>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> print(model)
For more examples, see the documentation of the methods below.
- Attributes:
- timestampstr
UTC timestamp at initialization.
- model_tdep
BaseModel
instance Time-dependent internal field model.
- model_static
BaseModel
instance Static internal field model.
- model_eulerdict of
Base
instances Dictionary containing the satellite’s name as key and the Euler angles as
Base
class instance.- model_caldict of
Base
instances Dictionary containing the satellite’s name as key and the calibration parameters as
Base
class instance.- n_smint, positive
Maximum spherical harmonic degree of external field in SM coordinates.
- coeffs_smndarray, shape (
n_sm
* (n_sm
+ 2),) Coefficients of static external field in SM coordinates.
- n_gsmint, positive
Maximum spherical harmonic degree of external field in GSM coordinates.
- coeffs_gsmndarray, shape (
n_gsm
* (n_gsm
+ 2),) Coefficients of static external field in GSM coordinates.
- breaks_deltadict with ndarrays, shape (\(m_q\) +1,)
Breaks of baseline corrections of static external field in SM coordinates. The dictionary keys are
'q10'
,'q11'
,'s11'
.- coeffs_deltadict with ndarrays, shape (1, \(m_q\))
Coefficients of baseline corrections of static external field in SM coordinates. The dictionary keys are
'q10'
,'q11'
,'s11'
.- namestr, optional
User defined name of the model.
- metadict, optional
Dictionary containing additional information about the model.
Methods
__call__
(time, radius, theta, phi[, rc_e, ...])Calculate the magnetic field of all sources from the CHAOS model.
from_mat
(filepath[, name, satellites])Alternative constructor for creating a
CHAOS
class instance.from_shc
(filepath, *[, name, leap_year])Alternative constructor for creating a
CHAOS
class instance.plot_maps_external
(time, radius, *[, nmax, ...])Plot global map of the external field from the CHAOS model.
plot_maps_static
(radius, *[, nmax])Plot global map of the static internal field from the CHAOS model.
plot_maps_tdep
(time, radius, *[, nmax, deriv])Plot global map of the time-dependent internal field from the CHAOS model.
plot_timeseries_tdep
(radius, theta, phi, ...)Plot the time series of the time-dependent internal field from the CHAOS model at a given location.
save_matfile
(filepath)Save CHAOS model to a mat-file.
save_shcfile
(filepath, *[, model, deriv, ...])Save spherical harmonic coefficients to a file in shc-format.
synth_coeffs_gsm
(time, *[, nmax, source])Compute the spherical harmonic coefficients of the far-magnetospheric magnetic field from the CHAOS model.
synth_coeffs_sm
(time, *[, nmax, source, rc])Compute the spherical harmonic coefficients of the near-magnetospheric magnetic field from the CHAOS model.
synth_coeffs_static
(*[, nmax])Compute the spherical harmonic coefficients of the static internal magnetic field from the CHAOS model.
synth_coeffs_tdep
(time, *[, nmax])Compute the spherical harmonic coefficients of the time-dependent internal magnetic field from the CHAOS model.
synth_euler_angles
(time, satellite, *[, ...])Extract Euler angles for specified satellite.
synth_values_gsm
(time, radius, theta, phi, *)Compute the vector components of the far-magnetospheric magnetic field from the CHAOS model.
synth_values_sm
(time, radius, theta, phi, *)Compute the vector components of the near-magnetospheric magnetic field from the CHAOS model.
synth_values_static
(radius, theta, phi, *[, ...])Compute the vector components of the static internal magnetic field from the CHAOS model.
synth_values_tdep
(time, radius, theta, phi, *)Compute the vector components of the time-dependent internal magnetic field from the CHAOS model.
- __call__(time, radius, theta, phi, rc_e=None, rc_i=None, source_list=None, nmax_static=None, verbose=None)[source]¶
Calculate the magnetic field of all sources from the CHAOS model.
All sources means the time-dependent and static internal fields, and the external magnetospheric (SM/GSM) fields including their induced parts.
- Parameters:
- timendarray, shape (…) or float
Array containing the time in modified Julian dates.
- radiusndarray, shape (…) or float
Radius of station in kilometers.
- thetandarray, shape (…) or float
Colatitude in degrees \([0^\circ, 180^\circ]\).
- phindarray, shape (…) or float
Longitude in degrees.
- rc_endarray, shape (…), optional
External part of the RC-index (defaults to linearly interpolating the hourly values given by the built-in RC-index file).
- rc_indarray, shape (…), optional
Internal part of the RC-index (defaults to linearly interpolating the hourly values given by the built-in RC-index file).
- source_listlist, [‘tdep’, ‘static’, ‘gsm’, ‘sm’] or str, {‘internal’, ‘external’}
Specify sources in any order. Default is all sources. Instead of a list, pass
source_list='internal'
which is equivalent tosource_list=['tdep', 'static']
(internal sources) orsource_list='external'
which is the same assource_list=['gsm', 'sm']
(external sources including induced part).- nmax_staticint, optional
Maximum spherical harmonic degree of the static internal magnetic field (defaults to 85).
- verbose{False, True}, optional
Print messages (defaults to
False
).
- Returns:
- B_radius, B_theta, B_phindarray, shape (…)
Radial, colatitude and azimuthal field components.
Examples
>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> Br, Bt, Bp = model(0., 6371.2, 45., 0., source_list=['tdep', 'static']) # only internal sources >>> Br array(-40418.23217586)
- classmethod from_mat(filepath, name=None, satellites=None)[source]¶
Alternative constructor for creating a
CHAOS
class instance.- Parameters:
- filepathstr
Path to mat-file containing the CHAOS model.
- namestr, optional
User defined name of the model. Defaults to the filename without the file extension.
- satelliteslist of strings, optional
List of satellite names whose Euler angles are stored in the mat-file. This is needed for correct referencing as this information is not given in the standard CHAOS mat-file format (defaults to
['oersted', 'champ', 'sac_c', 'swarm_a', 'swarm_b', 'swarm_c', 'cryosat-2_1', 'cryosat-2_2', 'cryosat-2_3']
.)
- Returns:
- model
CHAOS
instance Fully initialized model.
- model
See also
Examples
Load for example the mat-file
CHAOS-6-x7.mat
in the current working directory like this:>>> from chaosmagpy import CHAOS >>> model = CHAOS.from_mat('CHAOS-6-x7.mat') >>> print(model)
- classmethod from_shc(filepath, *, name=None, leap_year=None)[source]¶
Alternative constructor for creating a
CHAOS
class instance.- Parameters:
- filepathstr
Path to shc-file.
- namestr, optional
User defined name of the model. Defaults to
'CHAOS-<version>'
, where <version> is the default inbasicConfig['params.version']
.- 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.
- Returns:
- model
CHAOS
instance Class instance with either the time-dependent or static internal part initialized.
- model
See also
Examples
Load for example the shc-file
CHAOS-6-x7_tdep.shc
in the current working directory, containing the coefficients of time-dependent internal part of the CHAOS-6-x7 model.>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_shc('CHAOS-6-x7_tdep.shc') >>> print(model)
- plot_maps_external(time, radius, *, nmax=None, reference=None, source=None)[source]¶
Plot global map of the external field from the CHAOS model.
- Parameters:
- timendarray, shape (), (1,) or float
Time given as MJD2000 (modified Julian date).
- radiusndarray, shape (), (1,) or float
Radius in kilometers.
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- reference{‘all’, ‘GSM’, ‘SM’}, optional
Choose contribution either from GSM, SM or both added (default is ‘all’).
- source{‘all’, ‘external’, ‘internal’}, optional
Choose source to be external (inducing), internal (induced) or both added (default is ‘all’).
- Returns:
- B_radius, B_theta, B_phi
Global map of the radial, colatitude and azimuthal field components.
- plot_maps_static(radius, *, nmax=None, **kwargs)[source]¶
Plot global map of the static internal field from the CHAOS model.
- Parameters:
- radiusndarray, shape (), (1,) or float
Array containing the radius in kilometers.
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- **kwargskeywords
Other options are passed to
BaseModel.plot_maps()
method.
- Returns:
- B_radius, B_theta, B_phi
Global map of the radial, colatitude and azimuthal field components.
- plot_maps_tdep(time, radius, *, nmax=None, deriv=None, **kwargs)[source]¶
Plot global map of the time-dependent internal field from the CHAOS model.
- Parameters:
- timendarray, shape (), (1,) or float
Time given as MJD2000 (modified Julian date).
- radiusndarray, shape (), (1,) or float
Array containing the radius in kilometers.
- nmaxint, positive, optional
Maximum degree 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
BaseModel.plot_maps()
method.
- Returns:
- B_radius, B_theta, B_phi
Global map of the radial, colatitude and azimuthal field components.
- plot_timeseries_tdep(radius, theta, phi, **kwargs)[source]¶
Plot the time series of the time-dependent internal field from the CHAOS model at a given 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.
- nmaxint, positive, optional
Maximum degree 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.plot_timeseries()
method.
- Returns:
- B_radius, B_theta, B_phi
Time series plot of the radial, colatitude and azimuthal field components.
- save_matfile(filepath)[source]¶
Save CHAOS model to a mat-file.
The model must be fully specified as is the case if originally loaded from a mat-file.
- Parameters:
- filepathstr
Path and name of mat-file that is to be saved.
- save_shcfile(filepath, *, model=None, deriv=None, leap_year=None)[source]¶
Save spherical harmonic coefficients to a file in shc-format.
- Parameters:
- filepathstr
Path and name of output file *.shc.
- model{‘tdep’, ‘static’}, optional
Choose part of the model to save (default is ‘tdep’).
- derivint, optional
Derivative of the time-dependent field (default is 0, ignored for static source).
- 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.
- synth_coeffs_gsm(time, *, nmax=None, source=None)[source]¶
Compute the spherical harmonic coefficients of the far-magnetospheric magnetic field from the CHAOS model.
- Parameters:
- timendarray, shape (…) or float
Array containing the time in days.
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- source{‘external’, ‘internal’}, optional
Choose source either external or internal (default is ‘external’).
- Returns:
- coeffsndarray, shape (…,
nmax
* (nmax
+ 2)) Coefficients of the external GSM field in term of geographic coordinates (GEO).
- coeffsndarray, shape (…,
Examples
>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> model.synth_coeffs_gsm(0.0) array([11.63982782, -4.9276483 , -2.36281582, 0.46063709, -0.37934517, -0.18234297, 0.06281656, 0.07757099])
- synth_coeffs_sm(time, *, nmax=None, source=None, rc=None)[source]¶
Compute the spherical harmonic coefficients of the near-magnetospheric magnetic field from the CHAOS model.
- Parameters:
- timendarray, shape (…) or float
Array containing the time in days.
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- source{‘external’, ‘internal’}, optional
Choose source either external or internal (default is ‘external’).
- rcndarray, shape (…), optional
External (internal) part of the RC-index (defaults to linearly interpolating the hourly values given by the built-in RC-index file).
Note
Use the external part of the RC-index for
source='external'
(default) and the internal part forsource='internal'
.
- Returns:
- coeffsndarray, shape (…,
nmax
* (nmax
+ 2)) Coefficients of the external SM field coefficients in terms of geographic coordinates (GEO).
- coeffsndarray, shape (…,
Examples
>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> model.synth_coeffs_sm(0.) array([53.20309271, 3.79138724, -8.59458138, -0.62818711, 1.45506171, -0.57977672, -0.31660638, -0.43888236])
- synth_coeffs_static(*, nmax=None, **kwargs)[source]¶
Compute the spherical harmonic coefficients of the static internal magnetic field from the CHAOS model.
- Parameters:
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- **kwargskeywords
Other options are passed to
BaseModel.synth_coeffs()
method.
- Returns:
- coeffsndarray, shape (
nmax
* (nmax
+ 2),) Coefficients of the static internal field.
- coeffsndarray, shape (
Examples
>>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> model.synth_coeffs_static(nmax=50) array([ 0. , 0. , 0. , ..., 0.01655, -0.06339, 0.00715])
- synth_coeffs_tdep(time, *, nmax=None, **kwargs)[source]¶
Compute the spherical harmonic coefficients of the time-dependent internal magnetic field from the CHAOS model.
- Parameters:
- timendarray, shape (…) or float
Array containing the time in modified Julian dates.
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- **kwargskeywords
Other options to pass to
BaseModel.synth_coeffs()
method.
- Returns:
- coeffsndarray, shape (…,
nmax
* (nmax
+ 2)) Coefficients of the time-dependent internal field.
- coeffsndarray, shape (…,
Examples
>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> time = np.array([0., 10.])
>>> model.synth_coeffs_tdep(time, nmax=1) # dipole coefficients array([[-29614.72797782, -1728.47079907, 5185.50518939], [-29614.33800306, -1728.13680075, 5184.89196286]])
>>> model.synth_coeffs_tdep(time, nmax=1, deriv=1) # SV coefficients array([[ 14.25577646, 12.20214856, -22.43412895], [ 14.2317297 , 12.19625726, -22.36146885]])
- synth_euler_angles(time, satellite, *, dim=None, deriv=None, extrapolate=None)[source]¶
Extract Euler angles for specified satellite.
- Parameters:
- timefloat or ndarray, shape (…)
Time given as MJD2000 (modified Julian date).
- satellitestr
Satellite from which to get the euler angles.
- dimint, positive, optional
Maximum dimension (default is 3, for the three angles alpha, beta, gamma).
- derivint, positive, optional
Derivative in time (default is 0).
- extrapolate{‘linear’, ‘spline’, ‘constant’, ‘off’}, optional
Extrapolate to times outside of the model bounds. Defaults to
'linear'
.
- Returns:
- anglesndarray, shape (…, 3)
Euler angles alpha, beta and gamma in degrees, stored in trailing dimension.
Examples
>>> import chaosmagpy as cp >>> import numpy as np >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> time = np.array([500., 600.]) >>> model.meta['satellites'] # check satellite names ('oersted', 'champ', 'sac_c', 'swarm_a', 'swarm_b', 'swarm_c')
>>> model.synth_euler_angles(time, 'champ') array([[-0.05521985, -1.5763316 , 0.48787601], [-0.05440427, -1.57966925, 0.49043057]])
- synth_values_gsm(time, radius, theta, phi, *, nmax=None, source=None, grid=None)[source]¶
Compute the vector components of the far-magnetospheric magnetic field from the CHAOS model.
- Parameters:
- timefloat or ndarray, shape (…)
Time given as MJD2000 (modified Julian date).
- radiusfloat or ndarray, shape (…)
Array containing the radius in kilometers.
- thetafloat or ndarray, shape (…)
Array containing the colatitude in degrees \([0^\circ,180^\circ]\).
- phifloat or ndarray, shape (…)
Array containing the longitude in degrees.
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- source{‘all’, ‘external’, ‘internal’}, optional
Choose source to be external (inducing), internal (induced) or both added (default to ‘all’).
- gridbool, optional
If
True
, field components are computed on a regular grid, which is created fromtheta
andphi
as their outer product (defaults toFalse
).
- Returns:
- B_radius, B_theta, B_phindarray, shape (…)
Radial, colatitude and azimuthal field components.
Examples
>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> time = np.array([0., 10.]) >>> Br, Bt, Bp = model.synth_values_gsm(time, 6371.2, 45., 0.) >>> Br array([-8.18751916, -8.25661729])
- synth_values_sm(time, radius, theta, phi, *, rc_e=None, rc_i=None, nmax=None, source=None, grid=None)[source]¶
Compute the vector components of the near-magnetospheric magnetic field from the CHAOS model.
- Parameters:
- timefloat or ndarray, shape (…)
Time given as MJD2000 (modified Julian date).
- radiusfloat or ndarray, shape (…)
Array containing the radius in kilometers.
- thetafloat or ndarray, shape (…)
Array containing the colatitude in degrees \([0^\circ,180^\circ]\).
- phifloat or ndarray, shape (…)
Array containing the longitude in degrees.
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- source{‘all’, ‘external’, ‘internal’}, optional
Choose source to be external (inducing), internal (induced) or both added (default to ‘all’).
- gridbool, optional
If
True
, field components are computed on a regular grid, which is created fromtheta
andphi
as their outer product (defaults toFalse
).- rc_endarray, shape (…), optional
External part of the RC-index (defaults to linearly interpolating the hourly values given by the built-in RC-index file).
- rc_indarray, shape (…), optional
Internal part of the RC-index (defaults to linearly interpolating the hourly values given by the built-in RC-index file).
- Returns:
- B_radius, B_theta, B_phindarray, shape (…)
Radial, colatitude and azimuthal field components.
Examples
>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> time = np.array([0., 10.]) >>> Br, Bt, Bp = model.synth_values_sm(time, 6371.2, 45., 0.) >>> Br array([-18.85890361, -13.99893523])
- synth_values_static(radius, theta, phi, *, nmax=None, **kwargs)[source]¶
Compute the vector components of the static internal magnetic field from the CHAOS model.
- Parameters:
- radiusndarray, shape (…) or float
Radius of station in kilometers.
- thetandarray, shape (…) or float
Colatitude in degrees \([0^\circ, 180^\circ]\).
- phindarray, shape (…) or float
Longitude in degrees.
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- **kwargskeywords
Other options are passed to
BaseModel.synth_values()
method.
- Returns:
- B_radius, B_theta, B_phindarray, shape (…)
Radial, colatitude and azimuthal field components.
Examples
>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> Br, Bt, Bp = model.synth_values_static(6371.2, 45., 0., nmax=50) >>> Br array(-7.5608993)
- synth_values_tdep(time, radius, theta, phi, *, nmax=None, deriv=None, grid=None, extrapolate=None)[source]¶
Compute the vector components of the time-dependent internal magnetic field from the CHAOS model.
- Parameters:
- timendarray, shape (…) or float
Array containing the time in modified Julian dates.
- radiusndarray, shape (…) or float
Radius of station in kilometers.
- thetandarray, shape (…) or float
Colatitude in degrees \([0^\circ, 180^\circ]\).
- phindarray, shape (…) or float
Longitude in degrees.
- nmaxint, positive, optional
Maximum degree harmonic expansion (default is given by the model coefficients, but can also be smaller, if specified).
- derivint, positive, optional
Derivative in time (None defaults to 0). For secular variation, choose
deriv=1
.- gridbool, optional
If
True
, field components are computed on a regular grid. Arraystheta
andphi
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.
Examples
>>> import chaosmagpy as cp >>> model = cp.CHAOS.from_mat('CHAOS-6-x7.mat') >>> time = np.array([0., 10.])
Compute magnetic field components at specific location.
>>> Br, Bt, Bp = model.synth_values_tdep(time, 6371.2, 45., 0.) >>> Br array([-40422.44815265, -40423.15091334])
Only dipole contribution:
>>> Br, Bt, Bp = model.synth_values_tdep(time, 6371.2, 45., 0., nmax=1) >>> Br array([-44325.97679843, -44324.95294588])
Secular variation:
>>> Br, Bt, Bp = model.synth_values_tdep(time, 6371.2, 45., 0., deriv=1) >>> Br array([-25.64604374, -25.69002078])