chaosmagpy.model_utils.synth_values¶
- chaosmagpy.model_utils.synth_values(coeffs, radius, theta, phi, *, nmax=None, nmin=None, mmax=None, source=None, grid=None)[source]¶
Computes radial, colatitude and azimuthal field components from the magnetic potential field in terms of spherical harmonic coefficients.
- Parameters
- coeffsndarray, shape (…, M)
Coefficients of the spherical harmonic expansion. The last dimension is equal to the number of coefficients.
- 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 up to which expansion is to be used (default is given by the last dimension of
coeffs
, that is,M = nmax(nmax+2)
). This value can also be smaller if only coefficients at low degree should contribute.- nminint, positive, optional
Minimum degree of the expansion (defaults to 1). This value must correspond to the minimum degree of the coefficients and cannot be larger.
- mmaxint, positive, optional
Maximum order of the spherical harmonic expansion (defaults to
nmax
). This value can also be smaller. For example, if only the zonal part of the expansion should be used, setmmax = 0
.- source{‘internal’, ‘external’}, optional
Magnetic field source (default is an internal source).
- 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 (defaults toFalse
).
- Returns
- B_radius, B_theta, B_phindarray, shape (…)
Radial, colatitude and azimuthal field components.
Warning
The function can also evaluate the field components at the geographic poles, i.e. where
theta == 0.
ortheta == 180.
. However, users should be careful when doing this because the vector basis for spherical geocentric coordinates, \({{\mathbf{e}_r, \mathbf{e}_\theta, \mathbf{e}_\phi}}\), depends on longitude, which is not well defined at the poles. That is, at the poles, any value for the longitude maps to the same location in euclidean coordinates but gives a different vector basis in spherical geocentric coordinates. Nonetheless, by choosing a specific value for the longitude at the poles, users can define the vector basis, which then establishes the meaning of the spherical geocentric components returned by this function. The vector basis for the horizontal components is defined as\[\begin{split}\mathbf{e}_\theta &= \cos\theta\cos\phi\mathbf{e}_x - \cos\theta\sin\phi\mathbf{e}_y - \sin\theta\mathbf{e}_z\\ \mathbf{e}_\phi &= -\sin\phi\mathbf{e}_x + \cos\phi\mathbf{e}_y\end{split}\]Hence, at the geographic north pole as given by
theta = 0.
andphi = 0.
(chosen by the user), the returned componentB_theta
will be along the direction \(\mathbf{e}_\theta = \mathbf{e}_x\) andB_phi
along \(\mathbf{e}_\phi = \mathbf{e}_y\). However, ifphi = 180.
is chosen,B_theta
will be along \(\mathbf{e}_\theta = -\mathbf{e}_x\) andB_phi
along \(\mathbf{e}_\phi = -\mathbf{e}_y\).Notes
The function can work with different grid shapes, but the inputs have to satisfy NumPy’s broadcasting rules (
grid=False
, default). This also applies to the dimension of the coefficientscoeffs
excluding the last dimension.The optional parameter
grid
is for convenience. If set toTrue
, a singleton dimension is appended (prepended) totheta
(phi
) for broadcasting to a regular grid. The other inputsradius
andcoeffs
must then be broadcastable as before but now with the resulting regular grid.Examples
The most straight forward computation uses a fully specified grid. For example, compute the magnetic field at \(N=50\) grid points on the surface.
import chaosmagpy.model_utils as cpm import numpy as np N = 50 coeffs = np.ones((3,)) # degree 1 coefficients for all points radius = 6371.2 * np.ones((N,)) # radius of 50 points in km phi = np.linspace(-180., 180., num=N) # azimuth of 50 points in deg. theta = np.linspace(0., 180., num=N) # colatitude of 50 points in deg. B = cpm.synth_values(coeffs, radius, theta, phi) print([B[num].shape for num in range(3)]) # output shape (N,)
Instead of N points, compute the field on a regular \(N\times N\)-grid in azimuth and colatitude (slow).
radius_grid = 6371.2 * np.ones((N, N)) phi_grid, theta_grid = np.meshgrid(phi, theta) # grid of shape (N, N) B = cpm.synth_values(coeffs, radius_grid, theta_grid, phi_grid) print([B[num].shape for num in range(3)]) # output shape (N, N)
But this is slow since some computations on the grid are repeatedly done. The preferred method is to use NumPy’s broadcasting rules (fast).
radius_grid = 6371.2 # float, () or (1,)-shaped array broadcasted to NxN theta_grid = theta[:, None] # append singleton: (N, 1) phi_grid = phi[None, :] # prepend singleton: (1, N) B = cpm.synth_values(coeffs, radius_grid, theta_grid, phi_grid) print([B[num].shape for num in range(3)]) # output shape (N, N)
For convenience, you can do the same by using
grid=True
option, which appends a singleton dimension totheta
and inserts a singleton dimension at the second to last position in the shape ofphi
.B = cpm.synth_values(coeffs, radius_grid, theta, phi, grid=True) print([B[num].shape for num in range(3)]) # output shape (N, N)
Remember that
grid=False
(default) will result in (N,)-shaped outputs as in the first example.