pyrfu.pyrf package#

pyrfu.pyrf.autocorr(inp, maxlags: int | None = None, normed: bool = True)[source]#

Compute the autocorrelation function

Parameters:
  • inp (xarray.DataArray) – Input time series (scalar of vector).

  • maxlags (int, Optional) – Maximum lag in number of points. Default is None (i.e., len(inp) - 1).

  • normed (bool, Optional) – Flag to normalize the correlation.

Returns:

out – Autocorrelation function

Return type:

xarray.DataArray

pyrfu.pyrf.average_vdf(vdf, n_pts, method: str = 'mean')[source]#

Time averages the velocity distribution functions over n_pts in time.

Parameters:
  • vdf (xarray.DataArray) – Time series of the velocity distribution function.

  • n_pts (int) – Number of points (samples) of the averaging window.

  • method ({'mean', 'sum'}, Optional) – Method for averaging. Use ‘sum’ for counts. Default is ‘mean’.

Returns:

vdf_avg – Time series of the time averaged velocity distribution function.

Return type:

xarray.DataArray

pyrfu.pyrf.avg_4sc(b_list)[source]#

Computes the input quantity at the center of mass of the MMS tetrahedron.

Parameters:

b_list (list of xarray.DataArray) – List of the time series of the quantity for each spacecraft.

Returns:

b_avg – Time series of the input quantity a the enter of mass of the MMS tetrahedron.

Return type:

xarray.DataArray

Examples

>>> from pyrfu.mms import get_data
>>> from pyrfu.pyrf import avg_4sc

Time interval

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

Spacecraft indices

>>> b_mms = [get_data("B_gse_fgm_srvy_l2", tint, i) for i in range(1, 5)]
>>> b_xyz = avg_4sc(b_mms)
pyrfu.pyrf.c_4_grad(r_list, b_list, method: str = 'grad')[source]#

Calculate gradient of physical field using 4 spacecraft technique in [2] [3].

Parameters:
  • r_list (list of xarray.DataArray) – Time series of the positions of the spacecraft

  • b_list (list of xarray.DataArray) – Time series of the magnetic field at the corresponding positions

  • method ({"grad", "div", "curl", "bdivb", "curv"}, Optional) –

    Method flag :
    • ”grad” : compute gradient (default)

    • ”div” : compute divergence

    • ”curl” : compute curl

    • ”bdivb” : compute b.div(b)

    • ”curv” : compute curvature

Returns:

out – Time series of the derivative of the input field corresponding to the method

Return type:

xarray.DataArray

See also

pyrfu.pyrf.c_4_k

Calculates reciprocal vectors in barycentric coordinates.

References

[2]

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.

[3]

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. doi : https://www.issibern.ch/forads/sr-001-16.pdf

Examples

>>> from pyrfu.mms import get_data
>>> from pyrfu import mms, pyrf

Time interval

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

Load magnetic field and spacecraft position

>>> b_mms = [get_data("B_gse_fgm_srvy_l2", tint, m) for m in range(1, 5)]
>>> r_mms = [get_data("R_gse", tint, m) for m in range(1, 5)]
>>> gradb = pyrf.c_4_grad(r_mms, b_mms, "grad")
pyrfu.pyrf.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)
pyrfu.pyrf.c_4_k(r_list)[source]#

Calculates reciprocal vectors in barycentric coordinates.

Parameters:

r_list (list of xarray.DataArray) – Position of the spacecrafts.

Returns:

k_list – Reciprocal vectors in barycentric coordinates.

Return type:

list of xarray.DataArray

Notes

The units of reciprocal vectors are the same as [1/r].

pyrfu.pyrf.c_4_v(r_xyz, time)[source]#

Calculates velocity or time shift of discontinuity as in [6].

Parameters:
  • r_xyz (list) – Time series of the positions of the 4 spacecraft.

  • time (list) – Crossing times or time and velocity.

Returns:

out – Discontinuity velocity or time shift with respect to mms1.

Return type:

ndarray

References

[6]

Vogt, J., Haaland, S., and Paschmann, G. (2011) Accuracy of multi-point boundary crossing time analysis, Ann. Geophys., 29, 2239-2252, doi : https://doi.org/10.5194/angeo-29-2239-2011

pyrfu.pyrf.calc_ag(p_xyz)[source]#

Computes agyrotropy coefficient as in [16]

\[AG^{1/3} = \frac{|\operatorname[det]{\mathbf{P}} - \operatorname[det]{\mathbf{P}}|} {\operatorname[det]{\mathbf{P}} + \operatorname[det]{\mathbf{P}}}\]
Parameters:

p_xyz (xarray.DataArray) – Time series of the pressure tensor

Returns:

agyrotropy – Time series of the agyrotropy coefficient of the specie.

Return type:

xarray.DataArray

References

[16]

H. Che, C. Schiff, G. Le, J. C. Dorelli, B. L. Giles, and T. E. Moore (2018), Quantifying the effect of non-Larmor motion of electrons on the pres- sure tensor, Phys. Plasmas 25(3), 032101, doi: https://doi.org/10.1063/1.5016853.

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> ic = 1

Load magnetic field and electron pressure tensor

>>> b_xyz = mms.get_data("b_gse_fgm_srvy_l2", tint, 1)
>>> p_xyz_e = mms.get_data("pe_gse_fpi_fast_l2", tint, 1)

Rotate electron pressure tensor to field aligned coordinates

>>> p_fac_e_pp = mms.rotate_tensor(p_xyz_e, "fac", b_xyz, "pp")

Compute agyrotropy coefficient

>>> ag_e, ag_cr_e = pyrf.calc_ag(p_fac_e_pp)
pyrfu.pyrf.calc_agyro(p_xyz)[source]#

Computes agyrotropy coefficient as

\[A\Phi = \frac{|P_{\perp 1} - P_{\perp 2}|}{P_{\perp 1} + P_{\perp 2}}\]
Parameters:

p_xyz (xarray.DataArray) – Time series of the pressure tensor

Returns:

agyro – Time series of the agyrotropy coefficient of the specie.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> ic = 1

Load magnetic field and electron pressure tensor

>>> b_xyz = mms.get_data("b_gse_fgm_srvy_l2", tint, 1)
>>> p_xyz_e = mms.get_data("pe_gse_fpi_fast_l2", tint, 1)

Rotate electron pressure tensor to field aligned coordinates

>>> p_fac_e_qq = mms.rotate_tensor(p_xyz_e, "fac", b_xyz, "qq")

Compute agyrotropy coefficient

>>> agyro_e = pyrf.calc_agyro(p_fac_e_qq)
pyrfu.pyrf.calc_dng(p_xyz)[source]#

Computes agyrotropy coefficient as in [15]

\[D_{ng} = \frac{\sqrt{8 (P_{12}^2 + P_{13}^2 + P_{23}^2)}} {P_\parallel + 2 P_\perp}\]
Parameters:

p_xyz (xarray.DataArray) – Time series of the pressure tensor

Returns:

d_ng – Time series of the agyrotropy coefficient of the specie.

Return type:

xarray.DataArray

References

[15]

Aunai, N., M. Hesse, and M. Kuznetsova (2013), Electron nongyrotropy in the context of collisionless magnetic reconnection, Phys. Plasmas, 20(6), 092903, doi: https://doi.org/10.1063/1.4820953.

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> ic = 1

Load magnetic field and electron pressure tensor

>>> b_xyz = mms.get_data("b_gse_fgm_srvy_l2", tint, 1)
>>> p_xyz_e = mms.get_data("pe_gse_fpi_fast_l2", tint, 1)

Rotate electron pressure tensor to field aligned coordinates

>>> p_fac_e_pp = mms.rotate_tensor(p_xyz_e, "fac", b_xyz, "pp")

Compute agyrotropy coefficient

>>> d_ng_e = pyrf.calc_dng(p_fac_e_pp)
pyrfu.pyrf.calc_dt(inp)[source]#

Computes time step of the input time series.

Parameters:

inp (xarray.DataArray or xarray.Dataset) – Time series of the input variable.

Returns:

out – Time step in seconds.

Return type:

float

pyrfu.pyrf.calc_fs(inp)[source]#

Computes the sampling frequency of the input time series.

Parameters:

inp (xarray.DataArray or xarray.Dataset) – Time series of the input variable.

Returns:

out – Sampling frequency in Hz.

Return type:

float

pyrfu.pyrf.calc_sqrtq(p_xyz)[source]#

Computes agyrotropy coefficient as in [1]

\[Q = \frac{P_{12}^2 + P_{13}^2 + P_{23}^2} {P_\perp^2 + 2 P_\perp P_\parallel}\]
Parameters:

p_xyz (xarray.DataArray) – Time series of the pressure tensor

Returns:

sqrt_q – Time series of the agyrotropy coefficient of the specie

Return type:

xarray.DataArray

References

[1]

Swisdak, M. (2016), Quantifying gyrotropy in magnetic reconnection, Geophys. Res.Lett., 43, 43–49, doi: https://doi.org/10.1002/2015GL066980.

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> ic = 1

Load magnetic field and electron pressure tensor

>>> b_xyz = mms.get_data("b_gse_fgm_srvy_l2", tint, 1)
>>> p_xyz_e = mms.get_data("pe_gse_fpi_fast_l2", tint, 1)

Rotate electron pressure tensor to field aligned coordinates

>>> p_fac_e_pp = mms.rotate_tensor(p_xyz_e, "fac", b_xyz, "pp")

Compute agyrotropy coefficient

>>> sqrt_q_e = pyrf.calc_sqrtq(p_fac_e_pp)
pyrfu.pyrf.cart2sph(x, y, z)[source]#

Cartesian to spherical coordinate transform.

\[\begin{split}\alpha = \arctan \left( \frac{y}{x} \right) \\ \beta = \arccos \left( \frac{z}{r} \right) \\ r = \sqrt{x^2 + y^2 + z^2}\end{split}\]

with \(\alpha \in [0, 2\pi], \beta \in [0, \pi], r \geq 0\)

Parameters:
  • x (float or array_like) – x-component of Cartesian coordinates

  • y (float or array_like) – y-component of Cartesian coordinates

  • z (float or array_like) – z-component of Cartesian coordinates

Returns:

  • alpha (float or array_like) – Azimuth angle in radians

  • beta (float or array_like) – Elevation angle in radians (with 0 denoting North pole)

  • r (float or array_like) – Radius

pyrfu.pyrf.cart2sph_ts(inp, direction_flag: int = 1)[source]#

Computes magnitude, theta and phi angle from column vector xyz (first column is x ….) theta is 0 at equator. direction_flag = -1 -> to make transformation in opposite direction

Parameters:
  • inp (xarray.DataArray) – Time series to convert.

  • direction_flag ({1, -1}, Optional) – Set to 1 (default) to transform from cartesian to spherical coordinates. Set to -1 to transform from spherical to cartesian coordinates.

Returns:

out – Input field in spherical/cartesian coordinate system.

Return type:

xarray.DataArray

pyrfu.pyrf.cdfepoch2datetime64(epochs)[source]#

Converts CDF epochs to numpy.datetime64 with nanosecond precision.

Parameters:

epochs (float or int or array_like) – CDF epochs to convert.

Returns:

times – Array of times in datetime64([ns]).

Return type:

numpy.ndarray

pyrfu.pyrf.compress_cwt(cwt, nc: int = 100)[source]#

Compress the wavelet transform averaging of nc time steps.

Parameters:
  • cwt (xarray.Dataset) – Wavelet transform to compress.

  • nc (int, Optional) – Number of time steps for averaging. Default is 100.

Returns:

  • cwt_t (xarray.DataArray) – Sampling times.

  • cwt_x (ndarray) – Compressed wavelet transform of the first component of the field.

  • cwt_y (ndarray) – Compressed wavelet transform of the second component of the field.

  • cwt_z (ndarray) – Compressed wavelet transform of the third component of the field.

pyrfu.pyrf.convert_fac(inp, b_bgd, r_xyz: list | None = None)[source]#
Transforms to a field-aligned coordinate (FAC) system defined as :
  • R_parallel_z aligned with the background magnetic field

  • R_perp_y defined by R_parallel cross the position vector of the

    spacecraft (nominally eastward at the equator)

  • R_perp_x defined by R_perp_y cross R_par

If inp is one vector along r direction, out is inp[perp, para] projection.

Parameters:
Returns:

out – Time series of the input field in field aligned coordinates system.

Return type:

xarray.DataArray

Notes

All input parameters must be in the same coordinate system.

Examples

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

Time interval

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

Spacecraft index

>>> mms_id = 1

Load magnetic field (FGM) and electric field (EDP)

>>> b_xyz = mms.get_data("B_gse_fgm_brst_l2", tint, mms_id)
>>> e_xyz = mms.get_data("E_gse_edp_brst_l2", tint, mms_id)

Convert to field aligned coordinates

>>> e_xyzfac = pyrf.convert_fac(e_xyz, b_xyz, numpy.array([1, 0, 0]))
pyrfu.pyrf.corr_deriv(inp0, inp1, flag: bool = False)[source]#

Correlate the derivatives of two time series

Parameters:
  • inp0 (xarray.DataArray) – Time series of the first to variable to correlate with.

  • inp1 (xarray.DataArray) – Time series of the second to variable to correlate with.

  • flag (bool, Optional) – Flag if False (default) returns time instants of common highest first and second derivatives. If True returns time instants of common highest first derivative and zeros crossings.

Returns:

  • t1_d, t2_d (ndarray) – Time instants of common highest first derivatives.

  • t1_dd, t2_dd (ndarray) – Time instants of common highest second derivatives or zero crossings.

pyrfu.pyrf.cotrans(inp, flag, hapgood: bool = True)[source]#

Coordinate transformation GE0/GEI/GSE/GSM/SM/MAG as described in [1]

Parameters:
  • inp (xarray.DataArray or ndarray) – Time series of the input field.

  • flag (str) – Coordinates transformation “{coord1}>{coord2}”, where coord1 and coord2 can be geo/gei/gse/gsm/sm/mag.

  • hapgood (bool, Optional) – Indicator if original Hapgood sources should be used for angle computations or if updated USNO-AA sources should be used. Default = true, meaning original Hapgood sources.

Examples

>>> from pyrfu.mms import get_data
>>> from pyrfu.pyrf import cotrans

Time interval

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

Spacecraft index

>>> mms_id = 1

Load magnetic field in GSE coordinates

>>> b_gse = get_data("b_gse_fgm_srvy_l2", tint, mms_id)

Transform to GSM assuming that the original coordinates system is part of the inp metadata

>>> b_gsm = cotrans(b_gse, 'GSM')

If the original coordinates is not in the meta

>>> b_gsm = cotrans(b_gse, 'GSE>GSM')

Compute the dipole direction in GSE

>>> dipole = cotrans(b_gse.time, 'dipoledirectiongse')

References

[1]

Hapgood 1997 (corrected version of Hapgood 1992) Planet.Space Sci..Vol. 40, No. 5. pp. 71l - 717, 1992

pyrfu.pyrf.cross(inp1, inp2)[source]#

Computes cross product of two fields.

Parameters:
Returns:

out – Time series of the cross product Z = XxY.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Define time interval

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

Index of the MMS spacecraft

>>> mms_id = 1

Load magnetic field and electric field

>>> b_xyz = mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id)
>>> e_xyz = mms.get_data("E_gse_edp_fast_l2", tint, mms_id)

Compute magnitude of the magnetic field

>>> b_mag = pyrf.norm(b_xyz)

Compute ExB drift velocity

>>> v_xyz_exb = pyrf.cross(e_xyz, b_xyz) / b_mag ** 2
pyrfu.pyrf.date_str(tint, fmt: int = 1)[source]#

Creates a string corresponding to time interval for output plot naming.

Parameters:
  • tint (list of str) – Time interval.

  • fmt (int) –

    Format of the output :
    • 1 : “%Y%m%d_%H%M”

    • 2 : “%y%m%d%H%M%S”

    • 3 : “%Y%m%d_%H%M%S”_”%H%M%S”

    • 4 : “%Y%m%d_%H%M%S”_”%Y%m%d_%H%M%S”

Returns:

out – String corresponding to the time interval in the desired format.

Return type:

str

pyrfu.pyrf.datetime2iso8601(time)[source]#

Transforms datetime to TT2000 string format.

Parameters:

time (datetime.datetime) – Time to convert to tt2000 string.

Returns:

tt2000 – Time in TT20000 iso_8601 format.

Return type:

str

pyrfu.pyrf.datetime642iso8601(time)[source]#

Convert datetime64 in ns units to ISO 8601 time format .

Parameters:

time (ndarray) – Time in datetime64 in ns units.

Returns:

time_iso8601 – Time in ISO 8601 format.

Return type:

ndarray

pyrfu.pyrf.datetime642ttns(time)[source]#

Converts datetime64 in ns units to epoch_tt2000 (nanoseconds since J2000).

Parameters:

time (ndarray) – Times in datetime64 format.

Returns:

time_ttns – Times in epoch_tt2000 format (nanoseconds since J2000).

Return type:

ndarray

pyrfu.pyrf.datetime642unix(time)[source]#

Converts datetime64 in ns units to unix time.

Parameters:

time (ndarray) – Time in datetime64 format.

Returns:

time_unix – Time in unix format.

Return type:

ndarray

pyrfu.pyrf.dec_par_perp(inp, b_bgd, flag_spin_plane: bool = False)[source]#

Decomposes a vector into par/perp to B components. If flagspinplane decomposes components to the projection of b0 into the XY plane. alpha gives the angle between b0 and the XY. plane.

Parameters:
  • inp (xarray.DataArray) – Time series of the field to decompose.

  • b_bgd (xarray.DataArray) – Time series of the background magnetic field.

  • flag_spin_plane (bool, Optional) – Flag if True gives the projection in XY plane.

Returns:

  • a_para (xarray.DataArray) – Time series of the input field parallel to the background magnetic field.

  • a_perp (xarray.DataArray) – Time series of the input field perpendicular to the background magnetic field.

  • alpha (xarray.DataArray) – Time series of the angle between the background magnetic field and the XY plane.

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> mms_id = 1

Load magnetic field (FGM) and electric field (EDP)

>>> b_xyz = mms.get_data("B_gse_fgm_brst_l2", tint, mms_id)
>>> e_xyz = mms.get_data("E_gse_edp_brst_l2", tint, mms_id)

Decompose e_xyz into parallel and perpendicular to b_xyz components

>>> e_para, e_perp, _ = pyrf.dec_par_perp(e_xyz, b_xyz)
pyrfu.pyrf.dist_append(inp0, inp1)[source]#

Concatenate two distribution skymaps along the time axis.

Parameters:
Returns:

out – 3D skymap of the concatenated 3D skymaps.

Return type:

xarray.Dataset

Notes

The time series have to be in the correct time order.

pyrfu.pyrf.dot(inp1, inp2)[source]#

Computes dot product of two fields.

Parameters:
Returns:

out – Time series of the dot product Z = X.Y.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Load magnetic field, electric field and spacecraft position

>>> r_mms, b_mms, e_mms = [[] * 4 for _ in range(3)]
>>> for mms_id in range(1, 5):
>>>         r_mms.append(mms.get_data("R_gse", tint, mms_id))
>>>         b_mms.append(mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id))
>>>         e_mms.append(mms.get_data("E_gse_edp_fast_l2", tint, mms_id))

Compute current density using curlometer technique

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

Compute the electric at the center of mass of the tetrahedron

>>> e_xyz = pyrf.avg_4sc(e_mms)

Compute J.E dissipation

>>> je = pyrf.dot(j_xyz, e_xyz)
pyrfu.pyrf.dynamic_press(n_s, v_xyz, specie: str = 'ions')[source]#

Computes dynamic pressure.

Parameters:
  • n_s (xarray.DataArray) – Time series of the number density of the specie.

  • v_xyz (xarray.DataArray) – Time series of the bulk velocity of the specie.

  • specie ({"ions", "electrons"}, Optional) – Specie. Default “ions”.

Returns:

p_dyn – Time series of the dynamic pressure of the specie.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> mms_id = 1

Load ion bulk velocity and remove spintone

>>> v_xyz_i = mms.get_data("vi_gse_fpi_fast_l2", tint, mms_id)
>>> st_xyz_i = mms.get_data("sti_gse_fpi_fast_l2", tint, mms_id)
>>> v_xyz_i = v_xyz_i - st_xyz_i

Ion number density

>>> n_i = mms.get_data("ni_fpi_fast_l2", tint, mms_id)

Compute dynamic pressure

>>> p = pyrf.dynamic_press(n_i, v_xyz_i, specie="ions")
pyrfu.pyrf.e_vxb(v_xyz, b_xyz, flag: str = 'vxb')[source]#

Computes the convection electric field \(V\times B\) (default) or the \(E\times B/|B|^{2}\) drift velocity (flag=”exb”).

Parameters:
  • v_xyz (xarray.DataArray) – Time series of the velocity/electric field.

  • b_xyz (xarray.DataArray) – Time series of the magnetic field.

  • flag ({"vxb", "exb"}, Optional) –

    Method flag :
    • ”vxb” : computes convection electric field (Default).

    • ”exb” : computes ExB drift velocity.

Returns:

out – Time series of the convection electric field/ExB drift velocity.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> mms_id = 1

Load magnetic field and electric field

>>> b_xyz = mms.get_data("b_gse_fgm_srvy_l2", tint, mms_id)
>>> e_xyz = mms.get_data("e_gse_edp_fast_l2", tint, mms_id)

Compute ExB drift velocity

>>> v_xyz_exb = pyrf.e_vxb(e_xyz, b_xyz,"ExB")
pyrfu.pyrf.eb_nrf(e_xyz, b_xyz, v_xyz, flag=0)[source]#

Find E and B in MP system given B and MP normal vector.

Parameters:
  • e_xyz (xarray.DataArray) – Time series of the electric field.

  • b_xyz (xarray.DataArray) – Time series of the magnetic field.

  • v_xyz (xarray.DataArray) – Normal vector.

  • flag (str or ndarray) – Method flag : * a : L is along b_xyz, N closest to v_xyz and M = NxL * b : N is along v_xyz, L is the mean direction of b_xyz in plane perpendicular to N, and M = NxL * numpy,ndarray : N is along v_xyz , L is closest to the direction specified by L_vector (e.g., maximum variance direction), M = NxL

Returns:

out – to fill.

Return type:

xarray.DataArray

pyrfu.pyrf.ebsp(e_xyz, db_xyz, b_xyz, b_bgd, xyz, freq_int, **kwargs)[source]#

Calculates wavelet spectra of E&B and Poynting flux using wavelets (Morlet wavelet). Also computes polarization parameters of B using SVD [7]. SVD is performed on spectral matrices computed from the time series of B using wavelets and then averaged over a number of wave periods.

Parameters:
  • e_xyz (xarray.DataArray) – Time series of the wave electric field.

  • db_xyz (xarray.DataArray) – Time series of the wave magnetic field.

  • b_xyz (xarray.DataArray) – Time series of the high resolution background magnetic field used for E.B=0.

  • b_bgd (xarray.DataArray) – Time series of the background magnetic field used for field aligned coordinates.

  • xyz (xarray.DataArray) – Time series of the position time series of spacecraft used for field aligned coordinates.

  • freq_int (str or array_like) –

    Frequency interval :
    • ”pc12” : [0.1, 5.0]

    • ”pc35” : [2e-3, 0.1]

    • [fmin, fmax] : arbitrary interval [fmin,fmax]

  • polarization (bool) – Computes polarization parameters. Default False.

  • no_resample (bool) – No resampling, E and delta_b are given at the same time line. Default False.

  • fac (bool) – Uses FAC coordinate system (defined by b0 and optionally xyz), otherwise no coordinate system transformation is performed. Default True.

  • de_dot_b0 (bool) – Computes dEz from delta_b dot B = 0, uses full_b. Default False.

  • full_b_db (bool) – delta_b contains DC field. Default False.

  • nav (int) – Number of wave periods to average Default 8.

  • fac_matrix (numpy.ndarray) – Specify rotation matrix to FAC system Default None.

  • m_width_coeff (int or float) – Specify coefficient to multiple Morlet wavelet width by. Default 1.

Returns:

res

Dataset with :
  • txarray.DataArray

    Time.

  • fxarray.DataArray

    Frequencies.

  • bb_xxyyzzssxarray.DataArray
    delta_b power spectrum with :
    • […,0] : x

    • […,1] : y

    • […,2] : z

    • […,3] : sum

  • ee_xxyyzzssxarray.DataArray
    E power spectrum with :
    • […,0] : x

    • […,1] : y

    • […,2] : z

    • […,3] : sum

  • ee_ssxarray.DataArray

    E power spectrum (xx+yy spacecraft coordinates, e.g. ISR2).

  • pf_xyzxarray.DataArray

    Poynting flux (xyz).

  • pf_rtpxarray.DataArray

    Poynting flux (r, theta, phi) [angles in degrees].

  • dopxarray.DataArray

    3D degree of polarization.

  • dop2dxarray.DataArray

    2D degree of polarization in the polarization plane.

  • planarityxarray.DataArray

    Planarity of polarization.

  • ellipticityxarray.DataArray

    Ellipticity of polarization ellipse.

  • kxarray.DataArray

    k-vector (theta, phi FAC) [angles in degrees].

Return type:

xarray.Dataset

See also

pyrfu.plot.pl_ebsp

to fill.

pyrfu.pyrf.convert_fac

Transforms to a field-aligned coordinate.

Notes

This software was developed as part of the MAARBLE (Monitoring, Analyzing and Assessing Radiation Belt Energization and Loss) collaborative research project which has received funding from the European Community’s Seventh Framework Programme (FP7-SPACE-2011-1) under grant agreement n. 284520.

References

[7]

Santolík, O., Parrot. M., and Lefeuvre. F. (2003) Singular value decomposition methods for wave propagation analysis,Radio Sci., 38(1), 1010, doi : https://doi.org/10.1029/2000RS002523 .

Examples

>>> from pyrfu import mms, pyrf
>>> # Time interval
>>> tint_brst = ["2015-10-30T05:15:42.000", "2015-10-30T05:15:54.000"]
>>> # Spacecraft index
>>> mms_id = 3
>>> # Load spacecraft position
>>> tint_long = pyrf.extend_tint(tint_brst, [-100, 100])
>>> r_xyz = mms.get_data("R_gse", tint_long, mms_id)
>>> # Load background magnetic field, electric field and magnetic field
fluctuations
>>> b_xyz = mms.get_data("B_gse_fgm_brst_l2", tint_brst, mms_id)
>>> e_xyz = mms.get_data("E_gse_edp_brst_l2", tint_brst, mms_id)
>>> b_scm = mms.get_data("B_gse_scm_brst_l2", tint_brst, mms_id)
>>> # Polarization analysis
>>> options = dict(polarization=True, fac=True)
>>> polarization = pyrf.ebsp(e_xyz, b_scm, b_xyz, b_xyz, r_xyz,
>>>                          freq_int=[10, 4000], **options)
pyrfu.pyrf.edb(e_xyz, b_bgd, angle_lim: float = 20.0, flag_method: str = 'E.B=0')[source]#

Compute Ez under assumption \(\mathbf{E}.\mathbf{B}=0\) or \(\mathbf{E}.\mathbf{B} \approx 0\)

Parameters:
  • e_xyz (xarray.DataArray) – Time series of the electric field.

  • b_bgd (xarray.DataArray) – Time series of the background magnetic field.

  • angle_lim (float, Optional) – B angle with respect to the spin plane should be less than angle_lim degrees otherwise Ez is set to 0. Default is 20.

  • flag_method (str, Optional) –

    Assumption on the direction of the measured electric field :

    ”e.b=0” : \(\mathbf{E}.\mathbf{B}=0\). (Default) “e_par” : \(\mathbf{E}\) field along the B projection is coming from parallelelectric field. “e_perp+nan” : to fill.

Returns:

  • ed (xarray.DataArray) – Time series of the electric field output.

  • d (xarray.DataArray) – Time series of the B elevation angle above spin plane.

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft indices

>>> mms_id = 1

Load magnetic field, electric field and spacecraft position

>>> b_xyz = mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id)
>>> e_xyz = mms.get_data("E_gse_edp_fast_l2", tint, mms_id)

Compute Ez

>>> e_z, alpha = pyrf.edb(e_xyz, b_xyz)
pyrfu.pyrf.end(inp)[source]#

Gives the last time of the time series in unix format.

Parameters:

inp (xarray.DataArray or xarray.Dataset) – Time series of the input variable.

Returns:

out – Value of the last time in unix format.

Return type:

float

pyrfu.pyrf.estimate(what_to_estimate: str, radius: float, length: float | None = None)[source]#

Estimate values for some everyday stuff.

Parameters:
  • what_to_estimate (str) –

    Value to estimate:
    • ”capacitance_disk” estimates the capacitance of a disk

    (requires radius of the disk). * “capacitance_sphere” estimates of a sphere (requires radius of the sphere). * “capacitance_wire” estimates the capacitance of a wire (requires radius and length of the wire). * “capacitance_cylinder” estimates the capacitance of a cylinder (requires radius and half length of the cylinder).

  • radius (float) – Radius of the disk, sphere, wire or cylinder

  • length (float, Optional) – Length of the wire or half lenght of the cylinder.

Returns:

out – Estimated value.

Return type:

float

Examples

>>> from pyrfu import pyrf

Define radius of the sphere in SI units

>>> r_sphere = 20e-2

Computes the capacitance of the sphere

>>> c_sphere = pyrf.estimate("capacitance_sphere", r_sphere)
pyrfu.pyrf.extend_tint(tint, ext: list | None = None)[source]#

Extends time interval.

Parameters:
  • tint (list of str) – Reference time interval to extend.

  • ext (list of float or list of float) – Number of seconds to extend time interval [left extend, right extend].

Returns:

tint_new – Extended time interval.

Return type:

list of str

Examples

>>> from pyrfu import pyrf

Time interval

>>> tints = ["2015-10-30T05:15:42.000", "2015-10-30T05:15:54.000"]

Load spacecraft position

>>> tints_long = pyrf.extend_tint(tint, [-100, 100])
pyrfu.pyrf.filt(inp, f_min: float = 0.0, f_max: float = 1.0, order: int = -1)[source]#

Filters input quantity.

Parameters:
  • inp (xarray.DataArray) – Time series of the variable to filter.

  • f_min (float, Optional) – Lower limit of the frequency range. Default is 0. (Highpass filter).

  • f_max (float, Optional) – Upper limit of the frequency range. Default is 1. (Highpass filter).

  • order (int, Optional) – Order of the elliptic filter. Default is -1.

Returns:

out – Time series of the filtered signal.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

>>> tint = ["2017-07-18T13:03:34.000", "2017-07-18T13:07:00.000"]

Spacecraft index

>>> mms_id = 1

Load magnetic and electric fields

>>> b_xyz = mms.get_data("B_gse_fgm_brst_l2", tint, mms_id)
>>> e_xyz = mms.get_data("E_gse_edp_brst_l2", tint, mms_id)

Convert E to field aligned coordinates

>>> e_xyzfac = pyrf.convert_fac(e_xyz, b_xyz, [1,0,0])

Bandpass filter E waveform

>>> e_xyzfac_hf = pyrf.filt(e_xyzfac, 4, 0, 3)
>>> e_xyzfac_lf = pyrf.filt(e_xyzfac, 0, 4, 3)
pyrfu.pyrf.find_closest(inp1, inp2)[source]#

Finds pairs that are closest to each other in two time series.

Parameters:
  • inp1 (ndarray) – Vector with time instants.

  • inp2 (ndarray) – Vector with time instants.

Returns:

  • t1new (ndarray) – Identified time instants that are closest each other.

  • t2new (ndarray) – Identified time instants that are closest each other.

  • ind1new (ndarray) – Identified time instants that are closest each other.

  • ind2new (ndarray) – Identified time instants that are closest each other.

pyrfu.pyrf.get_omni_data(variables, tint, database: str = 'omni_hour')[source]#

Downloads OMNI data.

Parameters:
  • variables (list) – Keys of the variables to download.

  • tint (list) – Time interval.

  • database ({"omni_hour", "omni_min"}, Optional) – OMNI data resolution. Default is database = “omni_hour”.

Returns:

data – OMNI data.

Return type:

xarray.Dataset

pyrfu.pyrf.gradient(inp)[source]#

Computes time derivative of the input variable.

Parameters:

inp (xarray.DataArray) – Time series of the input variable.

Returns:

out – Time series of the time derivative of the input variable.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

>>> tint = ["2017-07-18T13:03:34.000", "2017-07-18T13:07:00.000"]

Spacecraft index

>>> mms_id = 1

Load magnetic field

>>> b_xyz = mms.get_data("B_gse_fgm_brst_l2", tint, mms_id)

Time derivative of the magnetic field

>>> db_dt = pyrf.gradient(b_xyz)
pyrfu.pyrf.gse2gsm(inp, flag: str = 'gse>gsm')[source]#

Converts GSE to GSM.

Parameters:
  • inp (xarray.DataArray or ndarray) – Time series of the input in GSE (GSM) coordinates. If ndarray first column is time in unix format.

  • flag ({"gse>gsm", "gsm>gse"}, Optional) – Flag for conversion direction. Default is “gse>gsm”

Returns:

out – Time series of the input in GSM (GSE) coordinates. If ndarray first column is time in unix format.

Return type:

xarray.DataArray or ndarray

See also

pyrfu.pyrf.geocentric_coordinate_transformation

pyrfu.pyrf.histogram(inp, bins=100, y_range=None, weights=None, density=True)[source]#

Computes 1D histogram of the inp with bins bins

Parameters:
  • inp (xarray.DataArray) – Time series of the input scalar variable.

  • bins (str or int or array_like, Optional) – Number of bins. Default is bins=100.

  • y_range ((float, float), Optional) – The lower and upper range of the bins. If not provided, range is simply (inp.min(), inp.max()). Values outside the range are ignored. The first element of the range must be less than or equal to the second. range affects the automatic bin computation as well. While bin width is computed to be optimal based on the actual data within range, the bin count will fill the entire range including portions containing no data.

  • weights (array_like, Optional) – An array of weights, of the same shape as inp. Each value in inp only contributes its associated weight towards the bin count (instead of 1). If density is True, the weights are normalized, so that the integral of the density over the range remains 1.

  • density (bool, Optional) – If False, the result will contain the number of samples in each bin. If True, the result is the value of the probability density function at the bin, normalized such that the integral over the range is 1. Note that the sum of the histogram values will not be equal to 1 unless bins of unity width are chosen; it is not a probability mass function.

Returns:

out – 1D distribution of the input time series.

Return type:

xarray.DataArray

pyrfu.pyrf.histogram2d(inp1, inp2, bins=100, y_range=None, weights=None, density=True)[source]#

Computes 2d histogram of inp2 vs inp1 with nbins number of bins.

Parameters:
  • inp1 (xarray.DataArray) – Time series of the x coordinates of the points to be histogrammed.

  • inp2 (xarray.DataArray) – Time series of the y coordinates of the points to be histogrammed.

  • bins (str or int or tuple, Optional) – Number of bins. Default is bins=100.

  • y_range (array_like, shape(2,2), Optional) – The leftmost and rightmost edges of the bins along each dimension (if not specified explicitly in the bins parameters): [[xmin, xmax], [ymin, ymax]]. All values outside of this range will be considered outliers and not tallied in the histogram.

  • weights (array_like, shape(N,), Optional) – An array of values w_i weighing each sample (x_i, y_i). Weights are normalized to 1 if normed is True. If normed is False, the values of the returned histogram are equal to the sum of the weights belonging to the samples falling into each bin.

  • density (bool, Optional) – If False, the default, returns the number of samples in each bin. If True, returns the probability density function at the bin, bin_count / sample_count / bin_area.

Returns:

out – 2D map of the density of inp2 vs inp1.

Return type:

xarray.DataArray

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_id = np.arange(1, 5)

Load magnetic field and electric field

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

Compute current density, etc

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

Compute magnitude of B and J

>>> b_mag = pyrf.norm(b_xyz)
>>> j_mag = pyrf.norm(j_xyz)

Histogram of J vs B

>>> h2d_b_j = pyrf.histogram2d(b_mag, j_mag)
pyrfu.pyrf.increments(inp, scale: int = 10)[source]#

Returns the increments of a time series.

\[y = |x_i - x_{i+s}|\]

where \(s\) is the scale.

Parameters:
  • inp (xarray.DataArray) – Input time series.

  • scale (int, Optional) – Scale at which to compute the increments. Default is 10.

Returns:

  • kurt (ndarray) – kurtosis of the increments, one per product, using the Fisher’s definition (0 value for a normal distribution).

  • result (xarray.DataArray) – An xarray containing the time series increments, one per product in the original time series.

pyrfu.pyrf.int_sph_dist(vdf, speed, phi, theta, speed_grid, **kwargs)[source]#

Integrate a spherical distribution function to a line/plane.

Parameters:
  • vdf (numpy.ndarray) – Phase-space density skymap.

  • speed (numpy.ndarray) – Velocity of the instrument bins,

  • phi (numpy.ndarray) – Azimuthal angle of the instrument bins.

  • theta (numpy.ndarray) – Elevation angle of the instrument bins.

  • speed_grid (numpy.ndarray) – Velocity grid for interpolation.

  • **kwargs – Keyw

pyrfu.pyrf.integrate(inp, time_step: float | None = None)[source]#

Integrate time series.

Parameters:
  • inp (xarray.DataArray) – Time series of the variable to integrate.

  • time_step (float, Optional) – Time steps threshold. All time_steps larger than 3*time_step are assumed data gaps, default is that time_step is the smallest value of all time_steps of the time series.

Returns:

out – Time series of the time integrated input.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

>>> tint = ["2015-12-14T01:17:40.200", "2015-12-14T01:17:41.500"]

Spacecraft index

>>> mms_id = 1

Load magnetic field and electric field

>>> b_xyz = mms.get_data("B_gse_fgm_brst_l2", tint, mms_id)
>>> e_xyz = mms.get_data("E_gse_edp_brst_l2", tint, mms_id)

Convert electric field to field aligned coordinates

>>> e_xyzfac = pyrf.convert_fac(e_xyz, b_xyz, [1, 0, 0])
pyrfu.pyrf.iplasma_calc(output: bool = False, verbose: bool = True)[source]#

Interactive function to calcute plasma paramters.

Parameters:
  • output (bool, Optional) – Flag to return dict with output. Default is False.

  • verbose (bool, Optional) – Flag to print the results function. Default is True.

Returns:

out – Hashtable with everything.

Return type:

dict

pyrfu.pyrf.iso86012datetime(time)[source]#

Converts ISO 8601 time to datetime

Parameters:

time (ndarray or list or str) – Time

Returns:

time_datetime – Time in datetime format.

Return type:

list

pyrfu.pyrf.iso86012datetime64(time)[source]#

Convert ISO8601 time format to datetime64 in ns units.

Parameters:

time (ndarray) – Time in ISO 8601 format

Returns:

time_datetime64 – Time in datetime64 in ns units.

Return type:

ndarray

pyrfu.pyrf.iso86012timevec(time)[source]#

Convert ISO 8601 time string into time vector.

Parameters:

time (ndarray or list or str) – Time in ISO 8601 format YYYY-MM-DDThh:mm:ss.mmmuuunnn.

Returns:

time_vec – Time vector.

Return type:

list

pyrfu.pyrf.iso86012unix(time)[source]#

Converts time in iso format to unix

Parameters:

time (str or array_like of str) – Time.

Returns:

out – Time in unix format.

Return type:

float or list of float

pyrfu.pyrf.l_shell(r_xyz)[source]#

Compute spacecraft position L Shell for a dipole magnetic field according to IRGF.

Parameters:

r_xyz (xarray.DataArray) – Time series of the spacecraft position. Must have a “COORDINATES_SYSTEM” attributes.

Returns:

out – Time series of the spacecraft position L-Shell.

Return type:

xarray.DataArray

pyrfu.pyrf.lowpass(inp, f_cut, fhz)[source]#

Filter the data through low or highpass filter with max frequency f_cut and subtract from the original.

Parameters:
  • inp (xarray.DataArray) – Time series of the input variable.

  • f_cut (float) – Cutoff frequency.

  • fhz (float) – Sampling frequency.

Returns:

out – Time series of the filter data.

Return type:

xarray.DataArray

pyrfu.pyrf.magnetosphere(model: str = 'mp_shue1998', tint: list | None = None)[source]#

Returns the location of magnetopause.

Parameters:
  • model (str) – Model to use. Implemented: ‘mp_shue1998’, ‘bs’. Default is ‘mp_shue1998’.

  • tint (list) – Time interval.

Returns:

  • x_ (ndarray) – X location of the magnetopause.

  • y_ (ndarray) – Y location of the magnetopause.

Examples

>>> from pyrfu.pyrf import magnetosphere
>>> x_mp, y_mp = magnetosphere("mp_shue1998", 10, -2)
pyrfu.pyrf.match_phibe_dir(b_xyz, e_xyz, angles: ndarray | None = None, f: float | None = None)[source]#

Get propagation direction by matching dBpar and “phi”. Tries different propagation directions and finds the direction perpendicular to the magnetic field that gives the best correlation between the electrostatic potential and the parallel wave magnetic field according to

\[\int E \textrm{d}t = \frac{B_0}{ne \mu_0} B_{wave}\]
Parameters:
  • b_xyz (xarray.DataArray) – Time series of the magnetic field (to be filtered if f is given).

  • e_xyz (xarray.DataArray) – Time series of the electric field (to be filtered if f is given).

  • angles (array_like, Optional) – The angles in degrees to try (1-180 default)

  • f (float, Optional) – Filter frequency.

Returns:

  • x (ndarray) – Normal direction (size: n_triesx3).

  • y (ndarray) – Propagation direction.

  • z (ndarray) – Magnetic field direction.

  • corr_vec (ndarray) – Correlation vector.

  • int_e_dt (ndarray) – Potential.

  • b_z (ndarray) – Wave magnetic field in parallel direction.

  • b_0 (ndarray) – Mean magnetic field.

  • de_k (ndarray) – Wave electric field in propagation direction.

  • de_n (ndarray) – Wave electric field in propagation normal direction.

  • e_k (ndarray) – Electric field in propagation direction.

  • e_n (ndarray) – Electric field in propagation normal direction.

pyrfu.pyrf.match_phibe_v(b_0, b_z, int_e_dt, n, v)[source]#

Get propagation velocity by matching dBpar and phi. Used together with irf_match_phibe_dir.m.Finds best match in amplitude given, B0, dB_par, phi, propagation direction implied, for specified n and v given as vectors.Returns a matrix of correlations and the two potentials that were correlated.

Parameters:
  • b_0 (array_like) – Average background magnetic field.

  • b_z (array_like) – Parallel wave magnetic field.

  • int_e_dt (array_like) – Potential.

  • n (array_like) – Vector of densities

  • v (array_like) – Vector of velocities.

Returns:

  • corr_mat (numpy.ndarray) – Correlation matrix(nn x nv).

  • phi_b (numpy.ndarray) – B0 * dB_par / n_e * e * mu0

  • phi_e (numpy.ndarray) – int(E) dt * v(dl=-vdt = > -dl = vdt)

pyrfu.pyrf.mean(inp, r_xyz, b_xyz, dipole_axis: DataArray | None = None)[source]#

Put inp into mean field coordinates defined by position vector r and magnetic field b if earth magnetic dipole axis z is given then uses another algorithm (good for auroral passages)

Parameters:
Returns:

out – Input field in mean field coordinates.

Return type:

xarray.DataArray

pyrfu.pyrf.mean_bins(inp0, inp1, bins: int = 10)[source]#

Computes mean of values of y corresponding to bins of x.

Parameters:
  • inp0 (xarray.DataArray) – Time series of the quantity of corresponding to the bins.

  • inp1 (xarray.DataArray) – Time series of the quantity to compute the binned mean.

  • bins (int, Optional) – Number of bins.

Returns:

out

Dataset with :
  • binsxarray.DataArray

    bin values of the x variable.

  • dataxarray.DataArray

    Mean values of y corresponding to each bin of x.

  • sigmaxarray.DataArray

    Standard deviation.

Return type:

xarray.Dataset

Examples

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

Time interval

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

Spacecraft indices

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

Load magnetic field and electric field

>>> r_mms, b_mms = [[] * 4 for _ in range(2)]
>>> for mms_id in range(1, 5):
>>>         r_mms.append(mms.get_data("R_gse", tint, mms_id))
>>>         b_mms.append(mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id))
>>>

Compute current density, etc

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

Compute magnitude of B and J

>>> b_mag = pyrf.norm(b_xyz)
>>> j_mag = pyrf.norm(j_xyz)

Mean value of J for 10 bins of B

>>> m_b_j = pyrf.mean_bins(b_mag, j_mag)
pyrfu.pyrf.mean_field(inp, deg)[source]#

Estimates mean field xm and wave field xw using polynomial fit of order deg for the number of columns larger than 3 assume that first column is time.

Parameters:
  • inp (array_like) – Input data.

  • deg (int) – Degree of the fitting polynomial

Returns:

  • inp_mean (numpy.ndarray) – Mean field.

  • inp_wave (numpy.ndarray) – Wave field

pyrfu.pyrf.medfilt(inp, n_pts: int = 11)[source]#

Applies a median filter over npts points to inp.

Parameters:
  • inp (xarray.DataArray) – Time series of the input variable.

  • n_pts (int, Optional) – Number of points of median filter.

Returns:

out – Time series of the median filtered input variable.

Return type:

xarray.DataArray

Examples

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

Time interval

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

Spacecraft indices

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

Load magnetic field and electric field

>>> r_mms, b_mms = [[] * 4 for _ in range(2)]
>>> for mms_id in range(1, 5):
>>>         r_mms.append(mms.get_data("R_gse", tint, mms_id))
>>>         b_mms.append(mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id))
>>>

Compute current density, etc

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

Get J sampling frequency

>>> fs = pyrf.calc_fs(j_xyz)

Median filter over 1s

>>> j_xyz = pyrf.medfilt(j_xyz,fs)
pyrfu.pyrf.median_bins(inp0, inp1, bins: int = 10)[source]#

Computes median of values of y corresponding to bins of x

Parameters:
  • inp0 (xarray.DataArray) – Time series of the quantity of bins.

  • inp1 (xarray.DataArray) – Time series of the quantity to the median.

  • bins (int, Optional) – Number of bins.

Returns:

out

Dataset with :
  • binsxarray.DataArray

    bin values of the x variable.

  • dataxarray.DataArray

    Median values of y corresponding to each bin of x.

  • sigmaxarray.DataArray

    Standard deviation.

Return type:

xarray.Dataset

Examples

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

Time interval

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

Spacecraft indices

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

Load magnetic field and electric field

>>> r_mms, b_mms = [[] * 4 for _ in range(2)]
>>> for mms_id in range(1, 5):
>>>         r_mms.append(mms.get_data("R_gse", tint, mms_id))
>>>         b_mms.append(mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id))
>>>

Compute current density, etc

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

Compute magnitude of B and J

>>> b_mag = pyrf.norm(b_xyz)
>>> j_mag = pyrf.norm(j_xyz)

Median value of J for 10 bins of B

>>> med_b_j = pyrf.mean_bins(b_mag, j_mag)
pyrfu.pyrf.movmean(inp, n_pts: int = 100)[source]#

Computes running average of the inp over npts points.

Parameters:
  • inp (xarray.DataArray) – Time series of the input variable.

  • n_pts (int, Optional) – Number of points to average over.

Returns:

out – Time series of the input variable averaged over npts points.

Return type:

xarray.DataArray

Notes

Works also with 3D skymap distribution.

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> mms_id = 1

Load ion pressure tensor

>>> p_xyz_i = mms.get_data("Pi_gse_fpi_brst_l2", tint, mms_id)

Running average the pressure tensor over 10s

>>> fs = pyrf.calc_fs(p_xyz_i)
>>> p_xyz_i = pyrf.movmean(p_xyz_i, int(10 * fs))
pyrfu.pyrf.mva(inp, flag: str = 'mvar')[source]#

Compute the minimum variance frame.

Parameters:
  • inp (xarray.DataArray) – Time series of the quantity to find minimum variance frame.

  • flag ({"mvar", "<bn>=0", "td"}, Optional) – Constrain. Default is “mvar”.

Returns:

  • out (xarray.DataArray) – Time series of the input quantity in LMN coordinates.

  • l (numpy.ndarray) – Eigenvalues l[0] > l[1] > l[2].

  • lmn (numpy.ndarray) – Eigenvectors LMN coordinates.

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> mms_id = 1

Load magnetic field

>>> b_xyz = mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id)

Compute MVA frame

>>> b_lmn, lamb, frame = pyrf.mva(b_xyz)
pyrfu.pyrf.mva_gui(inp)[source]#

GUI to interactively perform minimum variance analysis (MVA) on time series data by selecting the time interval to apply MVA on. The return of this function is a callback to the GUI object and class attributes like the minimum variance direction vector are accesable through this callback by the method get_minvar().

Parameters:

inp (xarray.DataArray) – Time series of the quantity to load into GUI and perform MVA on.

Returns:

Returns MvaGui object to access attributes. In order to keep GUI responsive and interactive, a reference to this object is needed.

Return type:

mva_callback

pyrfu.pyrf.new_xyz(inp, trans_mat)[source]#

Transform the input field to the new frame.

Parameters:
  • inp (xarray.DataArray) – Time series of the input field in the original coordinate system.

  • trans_mat (array_like) – Transformation matrix.

Returns:

out – Time series of the input in the new frame.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft indices

>>> mms_id = 1

Load magnetic field and electric field

>>> b_xyz = mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id)
>>> e_xyz = mms.get_data("E_gse_edp_fast_l2", tint, mms_id)

Compute MVA frame

>>> b_lmn, l, mva = pyrf.mva(b_xyz)

Move electric field to the MVA frame

>>> e_lmn = pyrf.new_xyz(e_xyz, mva)
pyrfu.pyrf.norm(inp)[source]#

Computes the magnitude of the input field.

Parameters:

inp (xarray.DataArray) – Time series of the input field.

Returns:

out – Time series of the magnitude of the input field.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> mms_id = 1

Load magnetic field

>>> b_xyz = mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id)

Compute magnitude of the magnetic field

>>> b_mag = pyrf.norm(b_xyz)
pyrfu.pyrf.normalize(inp)[source]#

Normalizes the input field.

Parameters:

inp (xarray.DataArray) – Time series of the input field.

Returns:

out – Time series of the normalized input field.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

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

Spacecraft index

>>> mms_id = 1

Load magnetic field

>>> b_xyz = mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id)

Compute the normalized magnetic field

>>> b = pyrf.normalize(b_xyz)
pyrfu.pyrf.optimize_nbins_1d(x, n_min: int = 1, n_max: int = 100)[source]#

Estimates the number of bins for 1d histogram that minimizes the risk function in [1] , obtained by direct decomposition of the MISE following the method described in [2] .

Parameters:
  • x (xarray.DataArray) – Input time series

  • n_min (int, Optional) – Minimum number of bins. Default is 1.

  • n_max (int, Optional) – Maximum number of bins. Default is 100.

Returns:

opt_n_x – Number of bins that minimizes the cost function.

Return type:

int

References

[1]

Rudemo, M. (1982) Empirical Choice of Histograms and Kernel Density Estimators. Scandinavian Journal of Statistics, 9, 65-78.

[2]

Shimazaki H. and Shinomoto S., A method for selecting the bin size of a time histogram Neural Computation (2007) Vol. 19(6), 1503-1527

pyrfu.pyrf.optimize_nbins_2d(x, y, n_min: list | None = None, n_max: list | None = None)[source]#

Estimates the number of bins for 2d histogram that minimizes the risk function in [1] , obtained by direct decomposition of the MISE following the method described in [2] .

Parameters:
  • x (xarray.DataArray) – Input time series of the first variable.

  • y (xarray.DataArray) – Input time series of the second variable.

  • n_min (array_like, Optional) – Minimum number of bins for each time series. Default is [1, 1]

  • n_max (array_like, Optional) – Maximum number of bins for each time series. Default is [100, 100]

Returns:

  • opt_n_x (int) – Number of bins of the first variable that minimizes the cost function.

  • opt_n_y (int) – Number of bins of the second variable that minimizes the cost function.

References

[1]

Rudemo, M. (1982) Empirical Choice of Histograms and Kernel Density Estimators. Scandinavian Journal of Statistics, 9, 65-78.

[2]

Shimazaki H. and Shinomoto S., A method for selecting the bin size of a time histogram Neural Computation (2007) Vol. 19(6), 1503-1527

pyrfu.pyrf.pid_4sc(r_mms, v_mms, p_mms, b_mms)[source]#

Compute Pi-D term using definition of [1] as :

\[Pi-D = - \Pi_{ij}D_{ij}\]

with \(\Pi_{ij}\) the deviatoric part of the pressure tensor :

\[ \begin{align}\begin{aligned}\Pi_{ij} = P_{ij} - p\delta_{ij}\\p = \frac{1}{3}P_{ii}\end{aligned}\end{align} \]

and \(D_{ij}\) the deviatoric part of the strain tensor :

\[ \begin{align}\begin{aligned}D_{ij} = \frac{1}{2}\left ( \partial_i u_j + \partial_j u_i \right ) - \frac{1}{3}\theta\delta_{ij}\\\theta = \nabla . u\end{aligned}\end{align} \]
Parameters:
Returns:

pid – Time series of the Pi-D.

Return type:

xarray.DataArray

References

[1]

Yang, Y., Matthaeus, W. H., Parashar, T. N., Wu, P., Wan, M., Shi, Y., et al. (2017). Energy transfer channels and turbulence cascade in Vlasov-Maxwell turbulence. Physical Review E, 95, 061201. doi : https://doi.org/10.1103/PhysRevE.95.061201

pyrfu.pyrf.plasma_beta(b_xyz, p_xyz)[source]#

Computes plasma beta at magnetic field sampling

\[eta =\]

rac{P_{th}}{P_b}

where : \(P_b = B^2 / 2 \mu_0\)

Parameters:
b_xyzxarray.DataArray

Time series of the magnetic field.

p_xyzxarray.DataArray

Time series of the pressure tensor.

Returns:
betaxarray.DataArray

Time series of the plasma beta at magnetic field sampling.

pyrfu.pyrf.plasma_calc(b_xyz, t_i, t_e, n_i, n_e)[source]#

Computes plasma parameters including characteristic length and time scales.

Parameters:
Returns:

out

Dataset of the plasma parameters :
  • timexarray.DataArray

    Time.

  • Wpexarray.DataArray

    Time series of the electron plasma frequency [rad.s^{-1}].

  • Fpexarray.DataArray

    Time series of the electron plasma frequency [Hz].

  • Wcexarray.DataArray

    Time series of the electron cyclotron frequency [rad.s^{-1}].

  • Fcexarray.DataArray

    Time series of the electron cyclotron frequency [Hz].

  • Wppxarray.DataArray

    Time series of the ion plasma frequency [rad.s^{-1}].

  • Fppxarray.DataArray

    Time series of the ion plasma frequency [Hz].

  • Fcpxarray.DataArray

    Time series of the ion cyclotron frequency [Hz].

  • Fuhxarray.DataArray

    Time series of the upper hybrid frequency [Hz].

  • Flhxarray.DataArray

    Time series of the lower hybrid frequency [Hz].

  • Vaxarray.DataArray

    Time series of the Alfvèn velocity (ions) [m.s^{-1}].

  • Vaexarray.DataArray

    Time series of the Alfvèn velocity (electrons) [m.s^{-1}].

  • Vtexarray.DataArray

    Time series of the electron thermal velocity [m.s^{-1}].

  • Vtpxarray.DataArray

    Time series of the electron thermal velocity [m.s^{-1}].

  • Vtsxarray.DataArray

    Time series of the sound speed [m.s^{-1}].

  • gamma_exarray.DataArray

    Time series of the electron Lorentz factor.

  • gamma_pxarray.DataArray

    Time series of the electron Lorentz factor.

  • Lexarray.DataArray

    Time series of the electron inertial length [m].

  • Lixarray.DataArray

    Time series of the electron inertial length [m].

  • Ldxarray.DataArray

    Time series of the Debye length [m].

  • Ndxarray.DataArray

    Time series of the number of electrons in the Debye sphere.

  • Roexarray.DataArray

    Time series of the electron Larmor radius [m].

  • Ropxarray.DataArray

    Time series of the ion Larmor radius [m].

  • Rosxarray.DataArray

    Time series of the length associated to the sound speed [m].

Return type:

xarray.Dataset

Examples

>>> from pyrfu import mms, pyrf

Time interval

>>> tint = ["2015-10-30T05:15:20.000", "2015-10-30T05:16:20.000"]

Spacecraft index

>>> mms_id = 1

Load magnetic field, ion/electron temperature and number density

>>> b_xyz = mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id)
>>> t_xyz_i = mms.get_data("Ti_gse_fpi_fast_l2", tint, mms_id)
>>> t_xyz_e = mms.get_data("Te_gse_fpi_fast_l2", tint, mms_id)
>>> n_i = mms.get_data("Ni_fpi_fast_l2", tint, mms_id)
>>> n_e = mms.get_data("Ne_fpi_fast_l2", tint, mms_id)

Compute scalar temperature

>>> t_xyzfac_i = mms.rotate_tensor(t_xyz_i, "fac", b_xyz, "pp")
>>> t_xyzfac_e = mms.rotate_tensor(t_xyz_e, "fac", b_xyz, "pp")
>>> t_i = pyrf.trace(t_xyzfac_i)
>>> t_e = pyrf.trace(t_xyzfac_e)

Compute plasma parameters

>>> plasma_params = pyrf.plasma_calc(b_xyz, t_i, t_e, n_i, n_e)
pyrfu.pyrf.poynting_flux(e_xyz, b_xyz, b_hat)[source]#

Estimates Poynting flux at electric field sampling as

\[S = \frac{E \times B}{\mu_0}\]

if b0 is given project the Poynting flux along b0

Parameters:
Returns:

  • s (xarray.DataArray) – Time series of the Pointing flux.

  • s_z (xarray.DataArray) – Time series of the projection of the Pointing flux (only if b0).

  • int_s (xarray.DataArray) – Time series of the time integral of the Pointing flux (if b0 integral along b0).

pyrfu.pyrf.pres_anis(p_fac, b_xyz)[source]#

Compute pressure anisotropy factor:

\[\mu_0 \frac{P_\parallel - P_\perp}{|\mathbf{B}|^2}\]
Parameters:
  • p_fac (xarray.DataArray) – Time series of the pressure tensor in field aligne coordinates.

  • b_xyz (xarray.DataArray) – Time series of the background magnetic field.

Returns:

p_anis – Time series of the pressure anisotropy.

Return type:

xarray.DataArray

See also

pyrfu.mms.rotate_tensor

Rotates pressure or temperature tensor into another coordinate system.

Examples

>>> from pyrfu import mms, pyrf

Time interval

>>> tint = ["2015-10-30T05:15:20.000", "2015-10-30T05:16:20.000"]

Spacecraft index

>>> mms_id = 1

Load magnetic field, ion/electron temperature and number density

>>> b_xyz = mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id)
>>> p_xyz_i = mms.get_data("Pi_gse_fpi_fast_l2", tint, mms_id)

Transform pressure tensor to field aligned coordinates >>> p_fac_i = mms.rotate_tensor(p_xyz_i, “fac”, b_xyz)

Compute pressure anistropy

>>> p_anis = pyrf.pres_anis(p_xyz_i, b_xyz)
pyrfu.pyrf.psd(inp, n_fft: int = 256, n_overlap: int = 128, window: str = 'hamming', d_flag: str = 'constant', scaling: str = 'density')[source]#

Estimate power spectral density using Welch’s method.

Welch’s method [11] computes an estimate of the power spectral density by dividing the data into overlapping segments, computing a modified periodogram for each segment and averaging the periodograms.

Parameters:
  • inp (xarray.DataArray) – Time series of measurement values.

  • n_fft (int, Optional) – Length of the FFT used, if a zero padded FFT is desired. Default to 256.

  • n_overlap (int, Optional) – Number of points to overlap between segments. Default to 128.

  • window (str, Optional) – Desired window to use. It is passed to get_window to generate the window values, which are DFT-even by default. See “get_window” or a list of windows and required parameters. Default Hanning

  • d_flag (str, Optional) – Specifies how to detrend each segment. It is passed as the “type” argument to the”detrend” function. Default to “constant”.

  • scaling (str, Optional) – Selects between computing the power spectral density (‘density’) where Pxx has units of V**2/Hz and computing the power spectrum (“spectrum”) where “Pxx” has units of V**2, if x is measured in V and “fs” is measured in Hz. Default to ‘density’

Returns:

out – Power spectral density or power spectrum of inp.

Return type:

xarray.DataArray

References

[11]

P. Welch, “The use of the fast Fourier transform for the estimation of power spectra: A method based on time averaging over short, modified periodograms”, IEEE Trans. Audio Electroacoust. vol. 15, pp. 70-73, 1967.

pyrfu.pyrf.pvi(inp, scale: int = 10)[source]#

Returns the PVI of a time series.

\[y = \frac{|x_i - x_{i+s}|^2}{<|x_i - x_{i+s}|^2>}\]

where \(s\) is the scale.

Parameters:
  • inp (xarray.DataArray) – Input time series.

  • scale (int, Optional) – Scale at which to compute the PVI. Default is 10.

Returns:

values – An xarray containing the pvi of the original time series.

Return type:

xarray.DataArray

pyrfu.pyrf.pvi_4sc(b_mms)[source]#

Compute the Partial Variance of Increments (PVI) using the definition in [1] as

\[PVI_{ij}(t) = \sqrt{\frac{|\Delta B_{ij}(t)|^2} {\langle|\Delta B_{ij}|^2\rangle}}\]

where \(\Delta B_{ij}(t) = B_i(t) - B_i(t)\) is the magnetic field increments, the average \(\langle . \rangle\) is taken over the whole interval, and \(i\), \(j\) = 1,2,3,4 is the MMS spacecraft number.

In addition, computes, the rotation of the magnetic field between two spacecraft, i.e., magnetic field shear angle, as :

\[\alpha_{ij}(t) = cos^{-1} \frac{B_i(t) . B_j(t)}{|B_i(t)| |B_j(t)|}\]
Parameters:

b_mms (list of xarray.DataArray) – List of the time series of the background magnetic field for the 4 spacecraft.

Returns:

  • pvi_ij (xarray.DataArray) – Time series of the Partial Variance of Increments for the 6 pairs of spacecraft.

  • alpha_ij (xarray.DataArray) – Time series of the magnetic field shear angle for the 6 pairs of spacecraft.

References

[1]

Chasapis, A., Retinó, A., Sahraoui, F., Vaivads, A., Khotyaintsev, Yu. V., Sundkvist, D., et al. (2015) Thin current sheets and associated electron heating in turbulent space plasma. Astrophys. J. Lett. 804:L1. doi: https://doi.org/10.1088/2041-8205/804/1/L1

pyrfu.pyrf.read_cdf(path, tint)[source]#

Reads CDF files.

Parameters:
  • path (str) – String of the filename in .cdf containing the L2 data

  • tint (list) – Time interval

Returns:

out_dict – Hash table with fields contained in the .cdf file.

Return type:

dict

pyrfu.pyrf.remove_repeated_points(inp)[source]#

Remove repeated elements in DataArray or structure data. Important when using defatt products. Must have a time variable.

Parameters:

inp (xarray.DataArray or dict) – Time series of the input variable.

Returns:

out – Time series of the cleaned input variable.

Return type:

xarray.DataArray or dict

pyrfu.pyrf.resample(inp, ref, method: str = '', f_s: float | None = None, window: int | None = None, thresh: float = 0)[source]#

Resample inp to the time line of ref. If sampling of X is more than two times higher than Y, we average X, otherwise we interpolate X.

Parameters:
  • inp (xarray.DataArray or xarray.Dataset) – Time series to resample.

  • ref (xarray.DataArray) – Reference time line.

  • method (str, Optional) – Method of interpolation “spline”, “linear” etc. (default “linear”) if method is given then interpolate independent of sampling.

  • f_s (float, Optional) – Sampling frequency of the Y signal, 1/window.

  • window (int or float or ndarray, Optional) – Length of the averaging window, 1/fsample.

  • thresh (float, Optional) – Points above STD*THRESH are disregarded for averaging

Returns:

out – Resampled input to the reference time line using the selected method.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

>>> tint = ["2015-10-30T05:15:20.000", "2015-10-30T05:16:20.000"]

Spacecraft index

>>> mms_id = 1

Load magnetic field and electric field

>>> b_xyz = mms.get_data("e_gse_fgm_srvy_l2", tint, mms_id)
>>> e_xyz = mms.get_data("e_gse_edp_fast_l2", tint, mms_id)

Resample magnetic field to electric field sampling

>>> b_xyz = pyrf.resample(b_xyz, e_xyz)
pyrfu.pyrf.shock_normal(spec, leq90: bool = True)[source]#

Calculates shock normals with different methods. Normal vectors are calculated by methods described in [1] and references therein.

The data can be averaged values or values from the time series in matrix format. If series is from time series all parameters are calculated from a random upstream and a random downstream point. This can help set errorbars on shock angle etc. The time series input must have the same size (up- and downstream can be different), so generally the user needs to resample the data first.

Parameters:
  • spec (dict) –

    Hash table with:
    • b_u : Upstream magnetic field (nT).

    • b_d : Downstream magnetic field.

    • v_u : Upstream plasma bulk velocity (km/s).

    • v_d : Downstream plasma bulk velocity.

    • n_u : Upstream number density (cm^-3).

    • n_d : Downstream number density.

    • r_xyz : Spacecraft position in time series format of 1x3 vector. Optional.

    • d2u : Down-to-up, is 1 or -1. Optional.

    • dt_f : Time duration of shock foot (s). Optional.

    • f_cp : Reflected ion gyrofrequency (Hz). Optional.

    • n : Number of Monte Carlo particles. Optional, default is 100.

  • leq90 (bool, Optional) – Force angles to be less than 90 (default). For leq90 = 0, angles can be between 0 and 180 deg. For time series input and quasi-perp shocks,leq90 = 0 is recommended.

Returns:

out

Hash table with:
  • n : Hash table containing normal vectors (n always points toward the

upstream region). From data:

  • mc : Magnetic coplanarity (10.14)

  • vc : Velocity coplanarity (10.18)

  • mx_1 : Mixed method 1 (10.15), [2]

  • mx_2 : Mixed method 2 (10.16), [2]

  • mx_3 : Mixed method 3 (10.17), [2]

From models (only if r_xyz is included in spec):
  • theta_bn : Angle between normal vector and b_u, same fields as n.

  • theta_vn : Angle between normal vector and v_u, same fields as n.

  • v_shHash table containing shock velocities:
    • gt : Using shock foot thickness (10.32). [8]

    • mf : Mass flux conservation (10.29).

    • sb : Using jump conditions (10.33). [9]

    • mo : Using shock foot thickness

  • infoHash table containing some more info:
    • msh : Magnetic shear angle.

    • vsh : Velocity shear angle.

    • cmat : Constraints matrix with normalized errors.

    • sig : Scaling factor to fit shock models to sc position. Calculated

    from (10.9-10.13) in [1]

Return type:

dict

References

[1] (1,2)

Schwartz, S. J. (1998), Shock and Discontinuity Normals, Mach Numbers, and Related Parameters, ISSI Scientific Reports Series, vol. 1, pp. 249–270.

[2] (1,2,3)

Abraham-Shrauner, B. (1972), “Determination of magnetohydrodynamic shock normals”, Journal of Geophysical Research, vol. 77, no. 4, p. 736. doi:10.1029/JA077i004p00736.

[3]

Farris, M. H., Petrinec, S. M., and Russell, C. T. (1991), The thickness of the magnetosheath: Constraints on the polytropic index”, Geophysical Research Letters, vol. 18, no. 10, pp. 1821–1824. doi:10.1029/91GL02090.

[4]

Slavin, J. A. and Holzer, R. E. (1981), Solar wind flow about the terrestrial planets, 1. Modeling bow shock position and shape, Journal of Geophysical Research, vol. 86, no. A13, pp. 11401–11418. doi:10.1029/JA086iA13p11401.

[5]

Peredo, M., Slavin, J. A., Mazur, E., and Curtis, S. A. (1995), Three-dimensional position and shape of the bow shock and their variation with Alfvénic, sonic, and magnetosonic Mach numbers and interplanetary magnetic field orientation, Journal of Geophysical Research, vol. 100, no. A5, pp. 7907–7916. doi:10.1029/94JA02545.

[6] (1,2)

Fairfield, D. H. (1971), Average and unusual locations of the Earth’s magnetopause and bow shock, Journal of Geophysical Research, vol. 76, no. 28, p. 6700, 1971. doi:10.1029/JA076i028p06700.

[7]

Formisano, V. (1979), Orientation and Shape of the Earth’s Bow Shock in Three Dimensions, Planetary and Space Science, vol. 27, no. 9, pp. 1151–1161. doi:10.1016/0032-0633(79)90135-1.

[8]

Gosling, J. T. and Thomsen, M. F. (1985), Specularly reflected ions, shock foot thicknesses, and shock velocity determination in space, Journal of Geophysical Research, vol. 90, no. A10, pp. 9893–9896. doi:10.1029/JA090iA10p09893.

[9]

Smith, E. J. and Burton, M. E. (1988), Shock analysis: Three useful new relations, Journal of Geophysical Research, vol. 93, no. A4, pp. 2730–2734. doi:10.1029/JA093iA04p02730.

pyrfu.pyrf.shock_parameters(spec)[source]#

Calculate shock related plasma parameters.

Parameters:

spec (dict) – Hash table with parameters input with fixed names. After the parameter name, a name of the region can be given, e.g. “u” and “d”. All parameters except b are optional.

Returns:

out – Hash table with derived plasma parameters from hash table spec with measured plasma parameters.

Return type:

dict

pyrfu.pyrf.solid_angle(inp0, inp1, inp2)[source]#

Calculates the solid angle of three vectors making up a triangle in a unit sphere with the sign taken into account.

Parameters:
  • inp0 (ndarray) – First vector.

  • inp1 (ndarray) – Second vector.

  • inp2 (ndarray) – Third vector.

Returns:

angle – Solid angle.

Return type:

float

pyrfu.pyrf.sph2cart(azimuth, elevation, r)[source]#

Transform spherical to cartesian coordinates

Parameters:
  • azimuth (float or ndarray) – Azimuthal angle (phi).

  • elevation (float or ndarray) – Elevation angle (theta)

  • r (float or ndarray) – Radius.

Returns:

  • x (float or ndarray) – Cartesian x-axis coordinates.

  • y (float or ndarray) – Cartesian y-axis coordinates.

  • z (float or ndarray) – Cartesian z-axis coordinates

pyrfu.pyrf.start(inp)[source]#

Gives the first time of the time series.

Parameters:

inp (xarray.DataArray) – Time series.

Returns:

out – Value of the first time in the desired format.

Return type:

float

pyrfu.pyrf.struct_func(inp, scales, order)[source]#

Returns the structure function of a time series

\[y= \frac{1}{N-s}\sum_{i=1}^{N-s}(x_i - x_{i+s})^o\]

where \(s\) is the scale, and \(o\) is the order.

Parameters:
  • inp (xarray.DataArray) – Input time series.

  • scales (array_like) – A list or an array containing the scales to calculate.

  • order (int) – Order of the exponential of the structure function.

  • ncut (int, Optional) – Number of standard deviation to cut (Kiyani et al., XXXX)

Returns:

values – An xarray containing the structure functions, one per product in the original time series. The index coordinate contains the scale value, and the attribute ‘order’ keeps a record on the order used for its calculation.

Return type:

xarray.DataArray

pyrfu.pyrf.t_eval(inp, times)[source]#

Evaluates the input time series at the target time.

Parameters:
  • inp (xarray.DataArray) – Time series if the input to evaluate.

  • times (ndarray) – Times at which the input will be evaluated.

Returns:

out – Time series of the input at times t.

Return type:

xarray.DataArray

pyrfu.pyrf.time_clip(inp, tint)[source]#

Time clip the input (if time interval is TSeries clip between start and stop).

Parameters:
Returns:

out – Time series of the time clipped input.

Return type:

xarray.DataArray

pyrfu.pyrf.timevec2iso8601(time)[source]#

Convert time vector into ISO 8601 format YYYY-MM-DDThh:mm:ss.mmmuuunnn.

Parameters:

time (ndarray) – Time vector

Returns:

time_iso8601 – Time in ISO 8601 format.

Return type:

ndarray

pyrfu.pyrf.trace(inp)[source]#

Computes trace of the time series of 2nd order tensors.

Parameters:

inp (xarray.DataArray) – Time series of the input 2nd order tensor.

Returns:

out – Time series of the trace of the input tensor.

Return type:

xarray.DataArray

Examples

>>> from pyrfu import mms, pyrf

Time interval

>>> tint = ["2015-10-30T05:15:20.000", "2015-10-30T05:16:20.000"]

Spacecraft index

>>> mms_id = 1

Load magnetic field and ion temperature

>>> b_xyz = mms.get_data("B_gse_fgm_srvy_l2", tint, mms_id)
>>> t_xyz_i = mms.get_data("Ti_gse_fpi_fast_l2", tint, mms_id)

Rotate to ion temperature tensor to field aligned coordinates

>>> t_xyzfac_i = mms.rotate_tensor(t_xyz_i, "fac", b_xyz, "pp")

Compute scalar temperature

>>> t_i = pyrf.trace(t_xyzfac_i)
pyrfu.pyrf.ts_append(inp0, inp1)[source]#

Concatenate two time series along the time axis.

Parameters:
Returns:

out – Concatenated time series.

Return type:

xarray.DataArray

Notes

The time series must be in the correct time order.

pyrfu.pyrf.ts_scalar(time, data, attrs: dict | None = None)[source]#

Create a time series containing a 0th order tensor

Parameters:
  • time (numpy.ndarray) – Array of times.

  • data (numpy.ndarray) – Data corresponding to the time list.

  • attrs (dict, Optional) – Attributes of the data list.

Returns:

out – 0th order tensor time series.

Return type:

xarray.DataArray

pyrfu.pyrf.ts_skymap(time, data, energy, phi, theta, **kwargs)[source]#

Creates a skymap of the distribution function.

Parameters:
  • time (array_like) – List of times.

  • data (array_like) – Values of the distribution function.

  • energy (array_like) – Energy levels.

  • phi (array_like) – Azimuthal angles.

  • theta (array_like) – Elevation angles.

  • **kwargs

    Hash table of keyword arguments with :
    • energy0array_like

      Energy table 0 (odd time indices).

    • energy1array_like

      Energy table 1 (even time indices).

    • esteptablearray_like

      Time series of the stepping table between energies (burst).

Returns:

out – Skymap of the distribution function.

Return type:

xarray.Dataset

pyrfu.pyrf.ts_spectr(time, ener, data, comp_name: str = 'energy', attrs: dict | None = None)[source]#

Create a time series containing a spectrum

Parameters:
  • time (numpy.ndarray) – Array of times.

  • ener (numpy.ndarray) – Y value of the spectrum (energies, frequencies, etc.)

  • data (numpy.ndarray) – Data of the spectrum.

  • attrs (dict, Optional) – Attributes of the data list.

Returns:

out – Time series of a spectrum

Return type:

xarray.DataArray

pyrfu.pyrf.ts_tensor_xyz(time, data, attrs: dict | None = None)[source]#

Create a time series containing a 2nd order tensor.

Parameters:
  • time (ndarray) – Array of times.

  • data (ndarray) – Data corresponding to the time list.

  • attrs (dict, Optional) – Attributes of the data list.

Returns:

out – 2nd order tensor time series.

Return type:

xarray.DataArray

pyrfu.pyrf.ts_time(time, attrs: dict | None = None)[source]#

Creates time line in DataArray.

Parameters:

time (ndarray) – Input time line.

Returns:

out – Time series of the time line.

Return type:

xarray.DataArray

pyrfu.pyrf.ts_vec_xyz(time, data, attrs: dict | None = None)[source]#

Create a time series containing a 1st order tensor.

Parameters:
  • time (ndarray) – Array of times.

  • data (ndarray) – Data corresponding to the time list.

  • attrs (dict, Optional) – Attributes of the data list.

Returns:

out – 1st order tensor time series.

Return type:

xarray.DataArray

pyrfu.pyrf.ttns2datetime64(time)[source]#

Convert time in epoch_tt2000 (nanosedconds since J2000) to datetime64 in ns units.

Parameters:

time (ndarray) – Time in epoch_tt2000 (nanoseconds since J2000) format.

Returns:

time_datetime64 – Time in datetime64 format in ns units.

Return type:

ndarray

pyrfu.pyrf.unix2datetime64(time)[source]#

Converts unix time to datetime64 in ns units.

Parameters:

time (ndarray) – Time in unix format.

Returns:

time_datetime64 – Time in datetime64 format.

Return type:

ndarray

pyrfu.pyrf.vht(e, b, no_ez: bool = False)[source]#

Estimate velocity of the De Hoffman-Teller frame from the velocity estimate the electric field eht=-vht x b

Parameters:
  • e (xarray.DataArray) – Time series of the electric field.

  • b (xarray.DataArray) – Time series of the magnetic field.

  • no_ez (boolean, Optional) – If True assumed no Ez. Default is False.

Returns:

  • vht (numpy.ndarray) – De Hoffman Teller frame velocity [km/s].

  • vht (xarray.DataArray) – Time series of the electric field in the De Hoffman frame.

  • dv_ht (numpy.ndarray) – Error of De Hoffman Teller frame.

pyrfu.pyrf.wave_fft(inp, window, frame_overlap: float = 10.0, frame_length: float = 20.0, f_sampling: float | None = None)[source]#

Short-Time Fourier Transform.

Parameters:
  • inp (xarray.DataArray) – Time series of the one dimension data.

  • window (str) – Window function such as rectwin, hamming (default).

  • frame_overlap (float, Optional) – Length of each frame overlaps in second.

  • frame_length (float, Optional) – Length of each frame in second.

  • f_sampling (float, Optional) – Sampling frequency.

Returns:

  • spectrogram (ndarray) – Spectrogram of x.

  • time (ndarray) – Value corresponds to the center of each frame (x-axis) in sec.

  • frequencies (ndarray) – Vector of frequencies (y-axis) in Hz.

pyrfu.pyrf.wavelet(inp, **kwargs)[source]#

Computes wavelet spectrogram based on fast FFT algorithm. :param inp: Input quantity. :type inp: xarray.DataArray :param **kwargs:

Hash table of keyword arguments with :
  • fsint or float

    Sampling frequency of the input time series.

  • flist or ndarray

    Vector [f_min f_max], calculate spectra between frequencies f_min and f_max.

  • nfint or float

    Number of frequency bins.

  • wavelet_widthint or float

    Width of the Morlet wavelet. Default 5.36.

  • linearfloat

    Linear spacing between frequencies of df.

  • return_powerbool

    Set to True to return the power, False for complex wavelet transform. Default True.

  • cut_edgebool

    Set to True to set points affected by edge effects to NaN, False to keep edge affect points. Default True

Returns:

out – Wavelet transform of the input.

Return type:

xarray.DataArray or xarray.Dataset

pyrfu.pyrf.wavepolarize_means(b_wave, b_bgd, min_psd: float = 1e-25, nop_fft: int = 256)[source]#

Analysis the polarization of magnetic wave using “means” method

Parameters:
  • b_wave (xarray.DataArray) – Time series of the magnetic field from Search Coil Magnetometer (SCM).

  • b_bgd (xarray.DataArray) – Time series of the magnetic field from Flux Gate Magnetometer (FGM).

  • min_psd (float, Optional) – Threshold for the analysis (e.g 1.0e-7). Below this value, the SVD analysis is meaningless if min_psd is not given, SVD analysis will be done for all waves. Default min_psd = 1e-25.

  • nop_fft (int, Optional) – Number of points in FFT. Default is 256.

Returns:

  • b_psd (xarray.DataArray) – Power spectrum density of magnetic filed wave.

  • wave_angle (xarray.DataArray) – (form 0 to 90)

  • deg_pol (xarray.DataArray) – Spectrogram of the degree of polarization (form 0 to 1).

  • elliptict (xarray.DataArray) – Spectrogram of the ellipticity (form -1 to 1)

  • helict (xarray.DataArray) – Spectrogram of the helicity (form -1 to 1)

Notes

b_wave and b_bgd should be from the same satellite and in the same coordinates

Warning

If one component is an order of magnitude or more greater than the other two then the polarization results saturate and erroneously indicate high degrees of polarization at all times and frequencies. Time series should be eyeballed before running the program. For time series containing very rapid changes or spikes the usual problems with Fourier analysis arise. Care should be taken in evaluating degree of polarization results. For meaningful results there should be significant wave power at the frequency where the polarization approaches 100%. Remember comparing two straight lines yields 100% polarization.

Examples

>>> from pyrfu import pyrf
>>> polarization = pyrf.wavepolarize_means(b_wave, b_bgd)
>>> polarization = pyrf.wavepolarize_means(b_wave, b_bgd, 1.0e-7)
>>> polarization = pyrf.wavepolarize_means(b_wave, b_bgd, 1.0e-7, 256)
pyrfu.pyrf.waverage(inp, f_sampl: float | None = None, n_pts: int = 7)[source]#

Computes weighted average.

Parameters:
  • inp (array_like) – Input data

  • f_sampl (float, Optional) – Sampling frequency.

  • n_pts (int, Optional) – Number of point to average over 5 ot 7. Default is 7

Returns:

out – Weighted averaged of inp

Return type:

ndarray

Submodules#