Multi-spacecraft analysis methods
This page demonstrates the use of multi-spacecraft analysis for MMS data. For more details about the package, see MultiSpacecraftAnalysis.jl.
using SPEDAS
using Dates
using CDAWeb
using DimensionalData
using CairoMakie, SpacePhysicsMakie
t0 = DateTime("2016-12-09T09:02")
t1 = DateTime("2016-12-09T09:04")
fgm_datasets = ntuple(4) do probe
get_dataset("MMS$(probe)_FGM_BRST_L2", t0, t1)
end
# Load data from CDF files into memory, MMS FGM data comes with magnitude
fields = ntuple(4) do probe
DimArray(fgm_datasets[probe]["mms$(probe)_fgm_b_gse_brst_l2"])[X(1:3), Ti(t0..t1)]
end
tplot(fields)
mec_datasets = ntuple(4) do probe
get_dataset("MMS$(probe)_MEC_BRST_L2_EPHT89D", t0, t1)
end
positions = ntuple(4) do probe
DimArray(mec_datasets[probe]["mms$(probe)_mec_r_gse"])[Ti(t0..t1)]
end
out = tlingradest(fields, positions)┌ 3×15358 DimStack ┐
├──────────────────┴───────────────────────────────────────────────────── dims ┐
↓ Y Sampled{Int64} Base.OneTo(3) ForwardOrdered Regular Points,
→ Ti Sampled{CommonDataFormat.TT2000} [2016-12-09T09:02:00.011, …, 2016-12-09T09:03:59.989] ForwardOrdered Irregular Points
├────────────────────────────────────────────────────────────────────── layers ┤
:Rbary eltype: Float64 dims: Y, Ti size: 3×15358
:Bbc eltype: Float64 dims: Y, Ti size: 3×15358
:Bmag eltype: Float64 dims: Ti size: 15358
:LGBx eltype: Float64 dims: Y, Ti size: 3×15358
:LGBy eltype: Float64 dims: Y, Ti size: 3×15358
:LGBz eltype: Float64 dims: Y, Ti size: 3×15358
:div eltype: Float64 dims: Ti size: 15358
:curl eltype: Float64 dims: Y, Ti size: 3×15358
:curv eltype: Float64 dims: Y, Ti size: 3×15358
:R_c eltype: Float64 dims: Ti size: 15358
└──────────────────────────────────────────────────────────────────────────────┘using SPEDAS: jparallel
jp = jparallel(out.Bbc, out.curl)
jp = setmeta(jp, ylabel = "Jparallel\n(nA/m²)")
tplot((out.Bbc, out.div, out.curl, out.curv, jp))
Validation with PySPEDAS
using PySPEDAS
using PythonCall
using Unitful
using Test
@py import pyspedas.projects.mms: mec, fgm, curlometer
trange = string.([t0, t1])
fgm_vars = @py fgm(probe=[1, 2, 3, 4], trange=trange, data_rate="brst", time_clip=true, varformat="*_gse_*")
mec_vars = @py mec(probe=[1, 2, 3, 4], trange=trange, data_rate="brst", time_clip=true, varformat="*_r_gse")
posits_py = ["mms1_mec_r_gse", "mms2_mec_r_gse", "mms3_mec_r_gse", "mms4_mec_r_gse"]
fields_py = ["mms1_fgm_b_gse_brst_l2", "mms2_fgm_b_gse_brst_l2", "mms3_fgm_b_gse_brst_l2", "mms4_fgm_b_gse_brst_l2"]
curlometer_vars = curlometer(fields=fields_py, positions=posits_py)
jp_py = PySPEDAS.get_data("jpar")
# Due to interpolation, pyspedas first and last values are NaN
@test (jp ./ u"A/m^2" .|> NoUnits) ≈ (jp_py[2:end-1]) atol=1e-5Test PassedBenchmark
using Chairmarks
@b tlingradest(fields, positions), curlometer(fields=fields_py, positions=posits_py)(81.241 ms (715 allocs: 7.210 MiB, 86.15% gc time), 5.251 s (69 allocs: 1.141 KiB))Julia is about 100 times faster than Python for similar workflows.
Dataset info
fgm_datasets[1]Dataset: ["/home/runner/.cdaweb/data/MMS1_FGM_BRST_L2/mms1_fgm_brst_l2_20161209090054_v5.87.0.cdf", "/home/runner/.cdaweb/data/MMS1_FGM_BRST_L2/mms1_fgm_brst_l2_20161209090304_v5.87.0.cdf"]
Group: mms1_fgm_brst_l2
Data variables
mms1_fgm_b_gse_brst_l2 (4 × 32000)
Datatype: Float32
Dimensions: nothing × Epoch
Attributes:
FILLVAL = Float32[-1.0f31]
FIELDNAM = Magnetic field vector in GSE plus Btotal (128 S/s)
VALIDMAX = Float32[20000.0, 20000.0, 20000.0, 20000.0]
SI_CONVERSION = 1.0e-9>T
TENSOR_ORDER = Int32[1]
CATDESC = Magnetic field vector in Geocentric Solar Ecliptic (GSE) cartesian coordinates plus Btotal (128 S/s)
DISPLAY_TYPE = time_series
VALIDMIN = Float32[-20000.0, -20000.0, -20000.0, 0.0]
UNITS = nT
REPRESENTATION_1 = represent_vec_tot
SCALETYP = linear
DEPEND_0 = Epoch
LABL_PTR_1 = StaticStrings.StaticString{10}["Bx GSE ", "By GSE ", "Bz GSE ", "Bt "]
COORDINATE_SYSTEM = GSE
FORMAT = E13.5
VAR_TYPE = data
mms1_fgm_b_gsm_brst_l2 (4 × 32000)
Datatype: Float32
Dimensions: nothing × Epoch
Attributes:
FILLVAL = Float32[-1.0f31]
FIELDNAM = Magnetic field vector in GSM plus Btotal (128 S/s)
VALIDMAX = Float32[20000.0, 20000.0, 20000.0, 20000.0]
SI_CONVERSION = 1.0e-9>T
TENSOR_ORDER = Int32[1]
CATDESC = Magnetic field vector in Geocentric Solar Magnetospheric (GSM) cartesian coordinates plus Btotal (128 S/s)
DISPLAY_TYPE = time_series
VALIDMIN = Float32[-20000.0, -20000.0, -20000.0, 0.0]
UNITS = nT
REPRESENTATION_1 = represent_vec_tot
SCALETYP = linear
DEPEND_0 = Epoch
LABL_PTR_1 = StaticStrings.StaticString{11}["Bx GSM ", "By GSM ", "Bz GSM ", "Bt "]
COORDINATE_SYSTEM = GSM
FORMAT = E13.5
VAR_TYPE = data
mms1_fgm_b_dmpa_brst_l2 (4 × 32000)
Datatype: Float32
Dimensions: nothing × Epoch
Attributes:
FILLVAL = Float32[-1.0f31]
FIELDNAM = Magnetic field vector in DMPA plus Btotal (128 S/s)
VALIDMAX = Float32[20000.0, 20000.0, 20000.0, 20000.0]
SI_CONVERSION = 1.0e-9>T
TENSOR_ORDER = Int32[1]
CATDESC = Magnetic field vector in Despun MPA-aligned cartesian coordinates plus Btotal (128 S/s)
DISPLAY_TYPE = time_series
VALIDMIN = Float32[-20000.0, -20000.0, -20000.0, 0.0]
UNITS = nT
VAR_NOTES = During nominal operatins in the region of interest, DMPA is within 3 degrees of GSE.
SCALETYP = linear
REPRESENTATION_1 = represent_vec_tot
DEPEND_0 = Epoch
LABL_PTR_1 = StaticStrings.StaticString{10}["Bx DMPA ", "By DMPA ", "Bz DMPA ", "Bt "]
COORDINATE_SYSTEM = DMPA>Despun Major Principal Axis
FORMAT = E13.5
VAR_TYPE = data
mms1_fgm_b_bcs_brst_l2 (4 × 32000)
Datatype: Float32
Dimensions: nothing × Epoch
Attributes:
FILLVAL = Float32[-1.0f31]
FIELDNAM = Magnetic field vector in BCS plus Btotal (128 S/s)
VALIDMAX = Float32[20000.0, 20000.0, 20000.0, 20000.0]
SI_CONVERSION = 1.0e-9>T
TENSOR_ORDER = Int32[1]
CATDESC = Magnetic field vector in Body Coordinate System cartesian coordinates plus Btotal (128 S/s)
DISPLAY_TYPE = time_series
VALIDMIN = Float32[-20000.0, -20000.0, -20000.0, 0.0]
UNITS = nT
REPRESENTATION_1 = represent_vec_tot
SCALETYP = linear
DEPEND_0 = Epoch
LABL_PTR_1 = StaticStrings.StaticString{10}["Bx BCS ", "By BCS ", "Bz BCS ", "Bt "]
COORDINATE_SYSTEM = BCS>Body Coordinate System
FORMAT = E13.5
VAR_TYPE = data
mms1_fgm_r_gse_brst_l2 (4 × 17)
Datatype: Float32
Dimensions: nothing × Epoch_state
Attributes:
FILLVAL = Float32[-1.0f31]
FIELDNAM = mms1_fgm_r_gse_brst_l2
VALIDMAX = Float32[1.0f6]
SI_CONVERSION = 1.0e3>m
TENSOR_ORDER = Int32[1]
CATDESC = Definitive Position in GSE coordinates, 30 second
DISPLAY_TYPE = time_series
VALIDMIN = Float32[-1.0f6]
UNITS = km
REPRESENTATION_1 = represent_vec_tot
SCALETYP = linear
DEPEND_0 = Epoch_state
LABL_PTR_1 = StaticStrings.StaticString{10}["X GSE ", "Y GSE ", "Z GSE ", "R "]
COORDINATE_SYSTEM = GSE
FORMAT = E12.2
VAR_TYPE = data
mms1_fgm_r_gsm_brst_l2 (4 × 17)
Datatype: Float32
Dimensions: nothing × Epoch_state
Attributes:
FILLVAL = Float32[-1.0f31]
FIELDNAM = mms1_fgm_r_gsm_brst_l2
VALIDMAX = Float32[1.0f6]
SI_CONVERSION = 1.0e3>m
TENSOR_ORDER = Int32[1]
CATDESC = Definitive Position in GSM coordinates, 30 second
DISPLAY_TYPE = time_series
VALIDMIN = Float32[-1.0f6]
UNITS = km
REPRESENTATION_1 = represent_vec_tot
SCALETYP = linear
DEPEND_0 = Epoch_state
LABL_PTR_1 = StaticStrings.StaticString{10}["X GSM ", "Y GSM ", "Z GSM ", "R "]
COORDINATE_SYSTEM = GSM
FORMAT = E12.2
VAR_TYPE = data
Support variables: Epoch, mms1_fgm_flag_brst_l2, Epoch_state, mms1_fgm_hirange_brst_l2, mms1_fgm_bdeltahalf_brst_l2, mms1_fgm_stemp_brst_l2, mms1_fgm_etemp_brst_l2, mms1_fgm_mode_brst_l2, mms1_fgm_rdeltahalf_brst_l2
Metadata variables: label_b_gse, label_b_gsm, label_b_dmpa, label_b_bcs, label_r_gse, label_r_gsm, represent_vec_tot
Global attributes
Project = ["STP>Solar-Terrestrial Physics"]
Source_name = ["MMS1>MMS Satellite Number 1"]
Discipline = ["Space Physics>Magnetospheric Science"]
Data_type = ["brst_l2"]
Descriptor = ["FGM>Flux Gate Magnetometer"]
File_naming_convention = ["source_descriptor_datatype_yyyyMMddHHmmss"]
Data_version = ["5.87.0"]
PI_name = ["J. Burch, C. Russell, W. Magnus"]
PI_affiliation = ["SWRI, UCLA, IWF"]
TEXT = ["The Fluxgate Magnetometers (FGM) on Magnetospheric Multiscale consist of a traditional Analog Fluxgate Magnetometer (AFG), and a Digital Fluxgate magnetometer (DFG). The dual magnetometers are operated as a single instrument providing a single intercalibrated data product. Range changes occur at different times on the two instruments so the gains checked each periapsis can be carried out unambiguously to apoapsis. Cross correlation of calibration parameters can separate causes of the any apparent calibration changes. Use of Electron Drift Instrument (EDI) to determine the field along the rotation axis allows accurate monitoring of the zero levels along the rotation axis. Prior to launch the magnetometers were calibrated at the Technical University, Braunschweig, except for the AFG magnetometers on MMS3 and MMS4, which were calibrated at UCLA. Both sets of sensors are operated for the entire MMS orbit, with slow survey (8 samples per second) outside of the Region of Interest (ROI), and fast survey (16 samples per second) inside the ROI. Within the ROI burst mode data (128 samples per second) are also acquired. A detailed description of the MMS fluxgate magnetometers, including science objectives, instrument description, calibration, magnetic cleanliness program, and data flow can be found at http://link.springer.com/article/10.1007%2Fs11214-014-0057-3 (DOI 10.1007/s11214-014-0057-3).Additional information can also be found at http://www-spc.igpp.ucla.edu/ssc/mms (UCLA),and http://www.iwf.oeaw.ac.at (IWF, Graz).", "For the purpose of creating a unified FGM Level2 data product, burst mode data is taken from DFG and survey mode data is taken from AFG. Because AFG and DFG are cross-calibrated on an orbit-averaged basis, small differences in offset may be observed between Level2 burst and survey mode data. Consequently, any differences are within the error of the measurement. Based on preliminary analysis of the data, the absolute error within the Region of Interest (ROI) is estimated to be no more than 0.1 nT in the spin-plane, 0.15 nT along the spin-axis and 0.2 nT in total magnitude."]
Instrument_type = ["Magnetic Fields (space)"]
Mission_group = ["MMS"]
Logical_source = ["mms1_fgm_brst_l2"]
Logical_file_id = ["mms1_fgm_brst_l2_20161209090054_v5.87.0"]
Logical_source_description = ["Level2 Flux Gate Magnetometer Burst DC Magnetic Field for MMS Satellite Number 1"]
Time_resolution = ["0.0078125 sec"]
Rules_of_use = ["See MMS Data Rights and Rules for Use https://lasp.colorado.edu/galaxy/display/mms/MMS+Data+Rights+and+Rules+for+Data+Use"]
Generated_by = ["UCLA"]
Generation_date = ["20170617"]
Acknowledgement = Any[]
MODS = ["version X=5: * Y-version number comes from cal file entries. \n * Ensures there are 2 ephemeris points before/after data to enable proper spline. \n * Fix to depend_0 of rdeltahalf: fixes bug when reading position data.\n * L-vector for DMPA2GSE transformation is smoothed with a gaussian filter, instead \n of using a single average value for the day. This short-term filter avoids \n introduding artificial jumps at 00:00 UTC and removes 7-minute 'wobble' after \n maneuvers in the GSE result. \n * Fixes error with DEFATT file selection found when choosing the \n daily DEFATT files to be used in Phase 2.\n * Fixed bug where reference Etemp was used for high range gain. Now uses measured Etemp.\nversion X=4: First version for public release of L2.\n Renamed variables to conform with new MMS variable name guidelines \n (obs_instr_paramName[_coordSys]_mode_level): \n Mag field parameters include 'b' for paramName. \n Use 'r' instead of 'pos' for S/C position paramName. \n Eliminated 'rate', replaced with 'bdeltahalf'. Added 'rdeltahalf'.\n l1a_mode is now just 'mode'.\nversion X=3: fixed removal of overlap between modes.\n fixed a bug that caused stemp and etemp to be empty.\nversion X=2: flag parameter name corrected: was 'status'\n added bits 4, 5, 6 to flag saturation on B1, B2, and B3, respectively\n added bit 7 to flag bad data at range changes\n Added etemp and l1a_mode parameters. \n rate, hirange, and stemp parameters now comply with MMS CDF Guidlelines, e.g.\n FILLVAL now defined for stemp and etemp, and is set to !values.f_nan\n No longer use Var_Parents attribute in stemp -- see Parents instead\n In this version, temperature-corrected gains are applied. Reference temperatures are used when \n stemp or etemp are set to FILLVAL. \n Non-linearity correction is applied to high rage DFG data.\nversion X=1: added 'flag', rate and hirange parameters (but 'flag' is actually called 'status')"]
ADID_ref = Any[]
LINK_TEXT = ["MMS home page", "SMART package home page", "Science Data Center", "FGM team home page", "FGM team"]
LINK_TITLE = ["at GSFC", "at SWRI", "at LASP", "at UCLA", "at IWF"]
HTTP_LINK = ["http://mms.gsfc.nasa.gov/", "http://mms.space.swri.edu/index.html", "https://lasp.colorado.edu/mms/sdc", "http://www-spc.igpp.ucla.edu/ssc/mms", "http://www.iwf.oeaw.ac.at"]
Parents = ["CDF>mms1_dfg_brst_l1a_20161209090054_v1.7.0", "CDF>mms1_fields_hk_l1b_10e_20161209_v0.5.0", "MMS1_DEFATT_2016343_2016344.V00", "MMS1_DEFEPH_2016343_2016344.V01", "CDF>mms1_dfg_lorangecal_l2pre_20100101_v12.87.0"]
References
- https://github.com/spedas/mms-examples/blob/master/basic/Curlometer%20Technique.ipynb