pyrfu.pyrf.c_4_j module#

pyrfu.pyrf.c_4_j.c_4_j(r_list, b_list)[source]#

Calculate current density \(J\) from using 4 spacecraft technique [1], the divergence of the magnetic field \(\nabla . B\), magnetic field at the center of mass of the tetrahedron, \(J \times B\) force, part of the divergence of stress associated with curvature \(\nabla.T_{shear}\) and gradient of the magnetic pressure \(\nabla P_b\). Where :

\[ \begin{align}\begin{aligned}J = \frac{\nabla \times B}{\mu_0}\\J \times B = \nabla.T_{shear} + \nabla P_b\\\nabla.T_{shear} = \frac{(B.\nabla) B}{\mu_0}\\\nabla P_b = \nabla \frac{B^2}{2\mu_0}\end{aligned}\end{align} \]

The divergence of the magnetic field is current density units as it shows the error on the estimation of the current density [2] .

Parameters:
Returns:

  • j (xarray.DataArray) – Time series of the current density [A.m^{-2}].

  • div_b (xarray.DataArray) – Time series of the divergence of the magnetic field [A.m^{-2}].

  • b_avg (xarray.DataArray) – Time series of the magnetic field at the center of mass of the tetrahedron, sampled at 1st SC time steps [nT].

  • jxb (xarray.DataArray) – Time series of the \(J\timesB\) force [T.A].

  • div_t_shear (xarray.DataArray) – Time series of the part of the divergence of stress associated with curvature units [T A/m^2].

  • div_pb (xarray.DataArray) – Time series of the gradient of the magnetic pressure.

See also

pyrfu.pyrf.c_4_k

Calculates reciprocal vectors in barycentric coordinates.

pyrfu.pyrf.c_4_grad

Calculate gradient of physical field using 4 spacecraft technique.

References

[1]

Dunlop, M. W., A. Balogh, K.-H. Glassmeier, and P. Robert (2002a), Four-point Cluster application of magnetic field analysis tools: The Curl- ometer, J. Geophys. Res., 107(A11), 1384, doi : https://doi.org/10.1029/2001JA005088.

[2]

Robert, P., et al. (1998), Accuracy of current determination, in Analysis Methods for Multi-Spacecraft Data, edited by G. Paschmann and P. W. Daly, pp. 395–418, Int. Space Sci. Inst., Bern. url : http://www.issibern.ch/forads/sr-001-16.pdf

Examples

>>> import numpy as np
>>> from pyrfu import mms, pyrf

Time interval

>>> tint = ["2019-09-14T07:54:00.000", "2019-09-14T08:11:00.000"]

Spacecraft indices

>>> mms_list = np.arange(1,5)

Load magnetic field and spacecraft position

>>> b_mms = [mms.get_data("b_gse_fgm_srvy_l2", tint, i) for i in mms_list]
>>> r_mms = [mms.get_data("r_gse_mec_srvy_l2", tint, i) for i in mms_list]

Compute current density, divergence of b, etc. using curlometer technique

>>> j_xyz, _, b_xyz, _, _, _ = pyrf.c_4_j(r_mms, b_mms)