Tutorial: Magnetic field processing¶

This tutorial will teach you how to perform the processing steps which require magnetic field calculations like mapping along the field line and calculations of adiabatic invariants.

We will use the time-binned variables from the cdf file for the tutorial. We also use energies and pitch angles for this tutorial and fold the pitch angles around 90 degrees assuming a symmetric pitch-angle distribution.

In [1]:

Copied!





import logging
import sys
from datetime import datetime, timedelta, timezone

from astropy import units as u

import el_paso as ep

ep.logger.setup_logging()


extraction_infos = [
    ep.ExtractionInfo(
        result_key="Epoch",
        name_or_column="Epoch_Ele",
        unit=ep.units.cdf_epoch,
    ),
    ep.ExtractionInfo(
        result_key="Energy",
        name_or_column="HOPE_ENERGY_Ele",
        unit=u.eV,
    ),
    ep.ExtractionInfo(
        result_key="Pitch_angle",
        name_or_column="PITCH_ANGLE",
        unit=u.deg,
        is_time_dependent=False,
    ),
    ep.ExtractionInfo(
        result_key="FEDU",
        name_or_column="FEDU",
        unit=(u.cm**2 * u.s * u.sr * u.keV) ** (-1),
    ),
    ep.ExtractionInfo(
        result_key="xGEO",
        name_or_column="Position_Ele",
        unit=u.km,
    ),
]

start_time = datetime(2017, 7, 14, tzinfo=timezone.utc)
end_time = datetime(2017, 7, 16, 23, 59, 59, tzinfo=timezone.utc)

file_name_stem = "rbspa_rel04_ect-hope-pa-l3_YYYYMMDD_.{6}.cdf"

ep.download(
    start_time,
    end_time,
    save_path=".",
    download_url="https://spdf.gsfc.nasa.gov/pub/data/rbsp/rbspa/l3/ect/hope/pitchangle/rel04/YYYY/",
    file_name_stem=file_name_stem,
    file_cadence="daily",
    method="request",
    skip_existing=True,
)

variables = ep.extract_variables_from_files(
    start_time, end_time, "daily", data_path=".", file_name_stem=file_name_stem, extraction_infos=extraction_infos
)

time_bin_methods = {
    "xGEO": ep.TimeBinMethod.NanMean,
    "Energy": ep.TimeBinMethod.NanMedian,
    "FEDU": ep.TimeBinMethod.NanMedian,
    "Pitch_angle": ep.TimeBinMethod.Repeat,
}

binned_time_variable = ep.processing.bin_by_time(
    variables["Epoch"],
    variables=variables,
    time_bin_method_dict=time_bin_methods,
    time_binning_cadence=timedelta(minutes=5),
)

variables["FEDU"].transpose_data([0, 2, 1])  # making it having dimensions (time, energy, pitch angle)
ep.processing.fold_pitch_angles_and_flux(variables["FEDU"], variables["Pitch_angle"])

# not needed anymore
del variables["Epoch"]
import logging
import sys
from datetime import datetime, timedelta, timezone

from astropy import units as u

import el_paso as ep

ep.logger.setup_logging()


extraction_infos = [
    ep.ExtractionInfo(
        result_key="Epoch",
        name_or_column="Epoch_Ele",
        unit=ep.units.cdf_epoch,
    ),
    ep.ExtractionInfo(
        result_key="Energy",
        name_or_column="HOPE_ENERGY_Ele",
        unit=u.eV,
    ),
    ep.ExtractionInfo(
        result_key="Pitch_angle",
        name_or_column="PITCH_ANGLE",
        unit=u.deg,
        is_time_dependent=False,
    ),
    ep.ExtractionInfo(
        result_key="FEDU",
        name_or_column="FEDU",
        unit=(u.cm**2 * u.s * u.sr * u.keV) ** (-1),
    ),
    ep.ExtractionInfo(
        result_key="xGEO",
        name_or_column="Position_Ele",
        unit=u.km,
    ),
]

start_time = datetime(2017, 7, 14, tzinfo=timezone.utc)
end_time = datetime(2017, 7, 16, 23, 59, 59, tzinfo=timezone.utc)

file_name_stem = "rbspa_rel04_ect-hope-pa-l3_YYYYMMDD_.{6}.cdf"

ep.download(
    start_time,
    end_time,
    save_path=".",
    download_url="https://spdf.gsfc.nasa.gov/pub/data/rbsp/rbspa/l3/ect/hope/pitchangle/rel04/YYYY/",
    file_name_stem=file_name_stem,
    file_cadence="daily",
    method="request",
    skip_existing=True,
)

variables = ep.extract_variables_from_files(
    start_time, end_time, "daily", data_path=".", file_name_stem=file_name_stem, extraction_infos=extraction_infos
)

time_bin_methods = {
    "xGEO": ep.TimeBinMethod.NanMean,
    "Energy": ep.TimeBinMethod.NanMedian,
    "FEDU": ep.TimeBinMethod.NanMedian,
    "Pitch_angle": ep.TimeBinMethod.Repeat,
}

binned_time_variable = ep.processing.bin_by_time(
    variables["Epoch"],
    variables=variables,
    time_bin_method_dict=time_bin_methods,
    time_binning_cadence=timedelta(minutes=5),
)

variables["FEDU"].transpose_data([0, 2, 1])  # making it having dimensions (time, energy, pitch angle)
ep.processing.fold_pitch_angles_and_flux(variables["FEDU"], variables["Pitch_angle"])

# not needed anymore
del variables["Epoch"]

[INFO    ] 2026-06-05 16:17:32 - el_paso.download:144 - Starting parallel download with 4 threads for 3 files...

[INFO    ] 2026-06-05 16:17:32 - el_paso.download:236 - File already exists, skipping download: rbspa_rel04_ect-hope-pa-l3_20170714_v7.4.0.cdf

[INFO    ] 2026-06-05 16:17:33 - el_paso.download:257 - Downloaded successfully: rbspa_rel04_ect-hope-pa-l3_20170715_v7.3.0.cdf

[INFO    ] 2026-06-05 16:17:33 - el_paso.download:257 - Downloaded successfully: rbspa_rel04_ect-hope-pa-l3_20170716_v7.3.0.cdf

[INFO    ] 2026-06-05 16:17:33 - el_paso.download:46 - download finished in 1.028 seconds

[INFO    ] 2026-06-05 16:17:33 - el_paso.extract_variables_from_files:112 - Extracting variables ...

[INFO    ] 2026-06-05 16:17:33 - el_paso.processing.bin_by_time:162 - Binning by time...

[WARNING ] 2026-06-05 16:17:34 - py.warnings:110 - /home/docs/checkouts/readthedocs.org/user_builds/el-paso/envs/latest/lib/python3.13/site-packages/el_paso/utils.py:209: UserWarning: Discarding nonzero nanoseconds in conversion.
  .to_pydatetime()

[WARNING ] 2026-06-05 16:17:35 - py.warnings:110 - /home/docs/checkouts/readthedocs.org/user_builds/el-paso/envs/latest/lib/python3.13/site-packages/el_paso/processing/bin_by_time.py:88: RuntimeWarning: All-NaN slice encountered
  binned_array = np.nanmedian(data, axis=0)

[INFO    ] 2026-06-05 16:17:35 - el_paso.processing.bin_by_time:68 - bin_by_time finished in 2.193 seconds

[INFO    ] 2026-06-05 16:17:35 - el_paso.processing.fold_pitch_angles_and_flux:98 - Folding pitch angles and flux ...

[WARNING ] 2026-06-05 16:17:35 - py.warnings:110 - /home/docs/checkouts/readthedocs.org/user_builds/el-paso/envs/latest/lib/python3.13/site-packages/el_paso/processing/fold_pitch_angles_and_flux.py:60: RuntimeWarning: Mean of empty slice
  folded_flux[:, :, i] = np.nanmean(masked_flux, axis=2)

[INFO    ] 2026-06-05 16:17:36 - el_paso.processing.fold_pitch_angles_and_flux:76 - fold_pitch_angles_and_flux finished in 0.049 seconds

Now we are ready to perform the magnetic field calculations. To make this most efficient, there exists a single function, which calculates all derived variables at once. The use specifies which variables should be calculated, which magnetic field model should be used, and which options for IRBEM are chosen (IRBEM is the underlying library performing the magnetic field calculations).

We will start by calculating just the variables which dependent on the satellite's position but not on the particles. The full list of available variable names can be found HERE.

In [2]:

Copied!





# Calculate magnetic field variables
from matplotlib import pyplot as plt

num_cores = 4
irbem_options = [1, 1, 4, 4, 0]

variables_to_compute: ep.processing.VariableRequest = [
    ("B_Calc", "T89"),
    ("MLT", "T89"),
    ("B_Eq", "T89"),
    ("R_Eq", "T89"),
]

magnetic_field_variables = ep.processing.compute_magnetic_field_variables(
    time_var=binned_time_variable,
    xgeo_var=variables["xGEO"],
    variables_to_compute=variables_to_compute,
    irbem_options=irbem_options,
    num_cores=num_cores,
)

plt.plot(binned_time_variable.get_data(), magnetic_field_variables["B_Calc_T89"].get_data(), "b", label="B_Calc")
plt.ylabel("B_local (nT)")

yyaxis = plt.gca().twinx()
yyaxis.plot(binned_time_variable.get_data(), magnetic_field_variables["R_Eq_T89"].get_data(), "r", label="R_Eq")
plt.xlabel("Time")
plt.ylabel("R_eq (RE)")
# Calculate magnetic field variables
from matplotlib import pyplot as plt

num_cores = 4
irbem_options = [1, 1, 4, 4, 0]

variables_to_compute: ep.processing.VariableRequest = [
    ("B_Calc", "T89"),
    ("MLT", "T89"),
    ("B_Eq", "T89"),
    ("R_Eq", "T89"),
]

magnetic_field_variables = ep.processing.compute_magnetic_field_variables(
    time_var=binned_time_variable,
    xgeo_var=variables["xGEO"],
    variables_to_compute=variables_to_compute,
    irbem_options=irbem_options,
    num_cores=num_cores,
)

plt.plot(binned_time_variable.get_data(), magnetic_field_variables["B_Calc_T89"].get_data(), "b", label="B_Calc")
plt.ylabel("B_local (nT)")

yyaxis = plt.gca().twinx()
yyaxis.plot(binned_time_variable.get_data(), magnetic_field_variables["R_Eq_T89"].get_data(), "r", label="R_Eq")
plt.xlabel("Time")
plt.ylabel("R_eq (RE)")

[INFO    ] 2026-06-05 16:17:36 - swvo.io.base:74 - KpOMNI data directory: /home/docs/.elpaso/OMNI_low_res

[INFO    ] 2026-06-05 16:17:36 - swvo.io.base:74 - KpNiemegk data directory: /home/docs/.elpaso/KpNiemegk

[INFO    ] 2026-06-05 16:17:36 - swvo.io.kp.read_kp_from_multiple_models:264 - Reading omni from 2017-07-14 00:00:22.425004+00:00 to 2017-07-16 23:55:22.425004+00:00

[INFO    ] 2026-06-05 16:17:36 - swvo.io.kp.read_kp_from_multiple_models:269 - Setting NaNs in omni from 2017-07-17 02:55:22.425004+00:00 to 2017-07-16 23:55:22.425004+00:00

[INFO    ] 2026-06-05 16:17:36 - swvo.io.utils:54 - Percentage of NaNs in data frame: 0.00%

[INFO    ] 2026-06-05 16:17:36 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:376 - 	Calculating local magnetic field values ...

[INFO    ] 2026-06-05 16:17:36 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:215 - get_local_B_field finished in 0.015 seconds

[INFO    ] 2026-06-05 16:17:36 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:320 - 	Calculating magnetic local time ...

[INFO    ] 2026-06-05 16:17:36 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:224 - get_MLT finished in 0.033 seconds

[INFO    ] 2026-06-05 16:17:36 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:133 - 	Calculating magnetic field and radial distance at the equator ...

  0%|          | 0/864 [00:00<?, ?it/s]

  0%|          | 0/864 [00:00<?, ?it/s]

  0%|          | 0/864 [00:01<?, ?it/s]

[INFO    ] 2026-06-05 16:17:37 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:227 - get_magequator finished in 1.090 seconds

Out[2]:

Text(0, 0.5, 'R_eq (RE)')

No description has been provided for this image

Next, we want to calculate equatorial pitch angles and adiabatic invariants. In this case, we have to provide the pitch-angles and energy variables as well.

The calculation of Lstar is quite heavy, so we will only calculate a small time frame for this tutorial.

In [3]:

Copied!





for var in variables.values():
    var.truncate(binned_time_variable, start_time=start_time, end_time=start_time + timedelta(hours=0.2))
binned_time_variable.truncate(binned_time_variable, start_time=start_time, end_time=start_time + timedelta(hours=0.2))

variables_to_compute: ep.processing.VariableRequest = [
    ("B_Calc", "T89"),
    ("MLT", "T89"),
    ("B_Eq", "T89"),
    ("R_Eq", "T89"),
    ("Alpha_Eq", "T89"),
    ("InvMu", "T89"),
    ("InvK", "T89"),
    ("L_star", "T89"),
]

magnetic_field_variables = ep.processing.compute_magnetic_field_variables(
    time_var=binned_time_variable,
    xgeo_var=variables["xGEO"],
    energy_var=variables["Energy"],
    pa_local_var=variables["Pitch_angle"],
    particle_species="electron",
    variables_to_compute=variables_to_compute,
    irbem_options=irbem_options,
    num_cores=1,
)

plt.plot(binned_time_variable.get_data(), magnetic_field_variables["R_Eq_T89"].get_data(), "b", label="R_eq")
plt.plot(binned_time_variable.get_data(), magnetic_field_variables["L_star_T89"].get_data(), "r", label="Lstar")

plt.xlabel("Time")
plt.legend()
plt.ylim((2, 6))
for var in variables.values():
    var.truncate(binned_time_variable, start_time=start_time, end_time=start_time + timedelta(hours=0.2))
binned_time_variable.truncate(binned_time_variable, start_time=start_time, end_time=start_time + timedelta(hours=0.2))

variables_to_compute: ep.processing.VariableRequest = [
    ("B_Calc", "T89"),
    ("MLT", "T89"),
    ("B_Eq", "T89"),
    ("R_Eq", "T89"),
    ("Alpha_Eq", "T89"),
    ("InvMu", "T89"),
    ("InvK", "T89"),
    ("L_star", "T89"),
]

magnetic_field_variables = ep.processing.compute_magnetic_field_variables(
    time_var=binned_time_variable,
    xgeo_var=variables["xGEO"],
    energy_var=variables["Energy"],
    pa_local_var=variables["Pitch_angle"],
    particle_species="electron",
    variables_to_compute=variables_to_compute,
    irbem_options=irbem_options,
    num_cores=1,
)

plt.plot(binned_time_variable.get_data(), magnetic_field_variables["R_Eq_T89"].get_data(), "b", label="R_eq")
plt.plot(binned_time_variable.get_data(), magnetic_field_variables["L_star_T89"].get_data(), "r", label="Lstar")

plt.xlabel("Time")
plt.legend()
plt.ylim((2, 6))

[INFO    ] 2026-06-05 16:17:37 - swvo.io.base:74 - KpOMNI data directory: /home/docs/.elpaso/OMNI_low_res

[INFO    ] 2026-06-05 16:17:37 - swvo.io.base:74 - KpNiemegk data directory: /home/docs/.elpaso/KpNiemegk

[INFO    ] 2026-06-05 16:17:37 - swvo.io.kp.read_kp_from_multiple_models:264 - Reading omni from 2017-07-14 00:00:22.425004+00:00 to 2017-07-14 00:10:22.425004+00:00

[INFO    ] 2026-06-05 16:17:37 - swvo.io.kp.read_kp_from_multiple_models:269 - Setting NaNs in omni from 2017-07-14 03:10:22.425004+00:00 to 2017-07-14 00:10:22.425004+00:00

[INFO    ] 2026-06-05 16:17:37 - swvo.io.utils:54 - Percentage of NaNs in data frame: 0.00%

[INFO    ] 2026-06-05 16:17:37 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:376 - 	Calculating local magnetic field values ...

[INFO    ] 2026-06-05 16:17:37 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:215 - get_local_B_field finished in 0.002 seconds

[INFO    ] 2026-06-05 16:17:37 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:320 - 	Calculating magnetic local time ...

[INFO    ] 2026-06-05 16:17:37 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:224 - get_MLT finished in 0.001 seconds

[INFO    ] 2026-06-05 16:17:37 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:133 - 	Calculating magnetic field and radial distance at the equator ...

  0%|          | 0/3 [00:00<?, ?it/s]

  0%|          | 0/3 [00:00<?, ?it/s]

  0%|          | 0/3 [00:01<?, ?it/s]

[INFO    ] 2026-06-05 16:17:38 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:227 - get_magequator finished in 1.030 seconds

[INFO    ] 2026-06-05 16:17:38 - el_paso.processing.compute_magnetic_field_variables:346 - 	Calculating equatorial pitch angle ...

[INFO    ] 2026-06-05 16:17:38 - el_paso.processing.compute_magnetic_field_variables:233 - Equatorial pitch angle calculation finished in 0.001 seconds

[INFO    ] 2026-06-05 16:17:38 - el_paso.processing.compute_magnetic_field_variables:398 - 	Calculating invariant mu ...

[INFO    ] 2026-06-05 16:17:38 - el_paso.processing.compute_magnetic_field_variables:242 - Invariant mu calculation finished in 0.001 seconds

[INFO    ] 2026-06-05 16:17:38 - el_paso.processing.compute_magnetic_field_variables:435 - 	Calculating invariant K ...

[INFO    ] 2026-06-05 16:17:38 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:588 - 	Calculating Lstar and J ...

  0%|          | 0/3 [00:00<?, ?it/s]

  0%|          | 0/3 [00:00<?, ?it/s]

  0%|          | 0/3 [00:01<?, ?it/s]

  0%|          | 0/3 [00:02<?, ?it/s]

 33%|███▎      | 1/3 [00:03<00:06,  3.00s/it]

 33%|███▎      | 1/3 [00:04<00:08,  4.01s/it]

 33%|███▎      | 1/3 [00:05<00:10,  5.01s/it]

 67%|██████▋   | 2/3 [00:06<00:03,  3.00s/it]

 67%|██████▋   | 2/3 [00:07<00:03,  3.50s/it]

 67%|██████▋   | 2/3 [00:08<00:04,  4.00s/it]

 67%|██████▋   | 2/3 [00:09<00:04,  4.51s/it]

[INFO    ] 2026-06-05 16:17:47 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:441 - get_Lstar finished in 9.036 seconds

[INFO    ] 2026-06-05 16:17:47 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:473 - 	Calculating mirror points ...

  0%|          | 0/3 [00:00<?, ?it/s]

  0%|          | 0/3 [00:00<?, ?it/s]

  0%|          | 0/3 [00:01<?, ?it/s]

[INFO    ] 2026-06-05 16:17:48 - el_paso.processing.magnetic_field_utils.magnetic_field_functions:443 - get_mirror_point finished in 1.027 seconds

[INFO    ] 2026-06-05 16:17:48 - el_paso.processing.compute_magnetic_field_variables:253 - Invariant K calculation finished in 10.066 seconds

Out[3]:

(2.0, 6.0)

In [ ]: