Statistical Characteristics of Energetic Electron Pitch Angle Distributions in the Van Allen Probe Era: 1. Butterfly Distributions With Flux Peaks at Preferred Pitch Angles

We present a statistical study on the properties of equatorial electron pitch angle distributions (PADs) in the Earth's outer radiation belt, for the first time based on particle measurements from the entire Van Allen Probes mission. A detailed selection criteria is used to identify intervals when flux measurements at energies from 0.2 to 3.4 MeV are available across a wide range of pitch angles close to the geomagnetic equatorial plane. To better characterize the shape of each pitch angle distribution, the flux data is fitted to functions of the equatorial pitch angle based on Legendre polynomials. Using this technique, we show that the shape of the PADs strongly depends on the particle's energy, location and geomagnetic activity. These results are used to identify the dominant physical processes responsible for creating PADs with different shapes. The results presented here mainly focus on the occurrence statistics and properties of butterfly PADs. Significantly, we find clear evidence that butterfly PADs have a peak flux that preferentially occurs at two distinct and discrete equatorial pitch angles, either 35° ± 5° or 65° ± 5°. Energy, L‐shell and magnetic local time variations indicate that the flux peaks at equatorial pitch angles ∼35° appear consistent with magnetopause shadowing, whilst those peaking at ∼65° may be associated with wave‐particle interactions. We also show that the flux at low and high pitch angles, for both butterfly and nonbutterfly PADs, remains remarkably well correlated even as the flux intensity varies by four orders of magnitude.


of 21
indicate that diffusive local acceleration of electrons from 0.1 to 3 MeV by magnetosonic waves would be most efficient for electrons with equatorial pitch angles between 50° and 80° and tend to produce PADs characterized by a local flux minimum at equatorial pitch angles near 90°, known as butterfly PADs (West et al., 1973). Ni et al. (2020) also presented observational evidence that at electron energies from 1.65 to 3.4 MeV butterfly PADs with increasing off-equatorial flux lasting >2 days are associated with the occurrence of magnetosonic waves. Omura et al. (2015) presented test particle simulations showing that electrons can be locally accelerated by a nonlinear interaction with chorus waves from keV to multi-MeV energies forming butterfly PADs with a peak flux at equatorial pitch angles between 50° and 80°. Xiao et al. (2015) showed that when the effects of diffusive local acceleration by both magnetosonic and chorus waves are combined, the formation of butterfly PADs observed in the outer radiation belt could be reproduced. Results presented in Albert et al. (2016) indicate that in the inner belt butterfly PADs can be formed by chorus and hiss waves when the effects of cross energy and pitch angle diffusion are included, without any magnetosonic waves present. Alternatively, the observed outer radiation belt flux enhancements may also result from adiabatic inward radial diffusion of the electrons by ultralow frequency (ULF) waves (see e,g., Fälthammar, 1966;Schulz & Lanzerotti, 1974). During inward radial diffusion, it is the equatorially mirroring electrons that are most efficiently accelerated by the ULF waves via a drift resonance (Southwood & Kivelson, 1982; see also Elkington et al., 1999;Fei et al., 2006;Ozeke et al., 2014), creating PADs peaked at equatorial pitch angles near 90°.
The loss of outer radiation belt electrons can be explained by either the particles entering the upper atmosphere or by the particles leaving the magnetosphere through the magnetopause, or by a combination of these two processes. Scattering of the particles' pitch angles by resonant interactions with a variety of waves can result in outer radiation belt flux loss to the upper atmosphere (see e.g., Reidy et al., 2021). Inside the plasmasphere (see Zhao et al., 2019; see also Breneman et al., 2015;Ripoll et al., 2019Ripoll et al., , 2021 it is seen that pitch angle scattering due to plasmaspheric hiss waves can occur, which strongly depends on the energy and pitch angle of the electrons, creating PADs that become highly peaked at 90° at lower electron energies, ∼100 keV. In addition to accelerating outer radiation belt electrons, chorus waves can also pitch angle scatter the particles into the loss cone, also changing the shape of the PADs ; see also Gu et al., 2012;Jin et al., 2018;Ma et al., 2018;Yang et al., 2016).
Theoretical studies suggest that a resonant interaction of the electrons with electromagnetic ion cyclotron (EMIC) waves can also rapidly scatter the pitch angles of energetic electrons into the loss cone causing them to enter the atmosphere (see e.g., Jordanova et al., 2008;Kersten et al., 2014;Ni et al., 2015;Ukhorskiy et al., 2010). Observational studies also suggest that EMIC waves may be responsible for the rapid storm time outer radiation belt flux dropouts observed on lower L-shells (see e.g., Xiang et al., 2018;Zhang et al., 2016). West et al. (1972) was the first to suggest that storm-time flux dropouts may result from loss to the magnetopause, via a process known as magnetopause shadowing. In addition, when the geomagnetic field is compressed on the dayside, drift-shell splitting causes the dayside particles with equatorial pitch angles closest to 90° to travel along asymmetric trajectories, which extend radially outward to greater distances on the dayside and reach the magnetopause first, creating butterfly PADs (Northrop & Teller, 1960;West, 1979;West et al., 1973; see also Sibeck et al., 1987). During storms with high levels of ULF wave power, magnetopause shadowing in combination with rapid outward radial diffusion can result in the loss of electrons at lower L-shells in the heart of the outer radiation belt (Turner et al., 2012; see also Loto'aniu et al., 2010;Mann & Ozeke, 2016;Ozeke et al., 2020).
Previous studies have analyzed the properties of electron PADs in the outer radiation (see e.g., Ni et al., 2016Ni et al., , 2020 as well as in the inner belt and slot regions (see e.g., Zhao et al., 2015; see also Albert, 2016;Zhao et al., 2014) to determine the possible physical processes responsible for creating the different PADs. In order to characterize the shape of each PAD, past studies have fitted the PAD to functions with the form sin n (α) where α is the pitch angle of the electrons and n is a constant (see Gannon et al., 2007; see also Allison et al., 2018;Shi et al., 2016). Butterfly PADs in general represent a small proportion of all observed PADs. However, determining the properties of these butterfly PADs and the conditions where they occur can reveal important information on the processes responsible for electron acceleration and loss in the outer radiation belt . However, these sin n (α) fits to the PADs are not able to represent the shape of the butterfly PADs. Chen et al. (2014) presented more detailed fits to the PADs based on Legendre polynomials, which are capable of better representing the shape of both butterfly and non-butterfly PADs (see also Zhao et al., 2021).

of 21
In this article, particle data from the entire Van Allen Probe mission in combination with the Legendre polynomial fits are used to investigate the properties of electron PADs in the outer radiation belt. The results presented here mainly focus on the occurrence statistics and properties of butterfly PADs, and the properties of nonbutterfly PADs will be investigated in more detail in future work. Here, we identify two distinct populations of butterfly PADs and investigate how the properties of these PADs depend on L-shell, magnetic local time (MLT), energy, and Dst, to help shed new light on the dominant physical processes responsible for creating these two different butterfly PAD populations.

Data Selection Criteria
A database of equatorial PADs is derived using electron flux and magnetic field measurements collected along the orbital paths of the NASA Van Allen Probes (Spence et al., 2013) from the start of the mission in September 2012 to the end of the mission in October 2019. Level 3 data from the Relativistic Electron-Proton Telescope (REPT) instrument  is used to determine the equatorial PADs at energies >2 MeV and level 3 flux data from the Magnetic Electron Ion Spectrometer (MagEIS) instrument  is used to determine the equatorial PADs at energies <2 MeV. Note that level 3 REPT data provides flux values in 17 local pitch angle bins, from 5.29° to 174.70° with a bin width of ±5.29°. Similarly, level 3 MagEIS data provides flux in 11 local pitch angle bins, from 8.18° to 171.81°, with a bin width of ±8.18. Combined REPT and MagEIS data products are also now available (see, Boyd et al., 2019Boyd et al., , 2021, although these new data products were not used in the current study. Equation 1 is used to derive the equatorial pitch angles of the particles from the local pitch angles measured by the Van Allen Probes, where it is assumed that the particles momentum and first adiabatic invariant are conserved as the particles move along the magnetic field line from the location of the probes to the geomagnetic equatorial plane (see e.g., Roederer & Zhang, 2016).
Here, B loc and α loc are the local magnetic field magnitude and particle pitch angle measured at the position of the probe. Similarly, B eq and α eq are the magnetic field magnitude and particle pitch angle along the same field line mapped to the geomagnetic equatorial plane. The TS04D magnetic field model (Tsyganenko & Sitnov, 2005) is used to determine both B eq and B loc .
In order to accurately characterize the shape of each PAD as a function of the equatorial pitch angle, data is only included in this study at times when the following conditions are satisfied: 1. Flux data is available at ≥9 equatorial pitch angles. 2. Flux data is available at ≥2 equatorial pitch angles within ±20° of 90°. 3. Flux data is available at equatorial pitch angles <35° and >145°T hese conditions ensure that there is flux data across a wide enough range of equatorial pitch angles to resolve the shape of each PAD. Note that the low 10.2⁰ inclination of the Van Allen Probes orbit ensures that only particle data close to the magnetic equatorial plane is used in this study. In addition, the second data selection condition listed above further restricts the particle data to locations close to the magnetic equatorial plane. Each individual equatorial PAD is binned with L*, at L* = 3.0, 3.5, 4.0, 4.5, and 5.0, with a width of ±0.25, derived using the TS04D magnetic field model. For simplicity all L*-values are derived at an equatorial pitch angle of 90°. Note that L*-values derived at high equatorial pitch angles can be significantly larger than those derived at low equatorial pitch angles (see e.g., Tu et al., 2019). However, in general, below L*∼5, L* is not a strong function of the equatorial pitch angle, as illustrated in Supplementary Figure S1 in Supporting Information S1. The equatorial PADs are also binned with Dst, at −60, −30, and 0 nT, with a bin width of ±15 nT.

Fitting the PADs
Each PAD, determined at the discrete equatorial pitch angles, is fit to a simple function of the equatorial pitch angle. We follow a similar approach to that outlined in Chen et al. (2014) (see also Zhao et al., 2021) and use a fitting function based on Legendre polynomials given by Here, log[j(α eq )] is the natural log of the differential electron flux as a function of the equatorial pitch angle, C n are the Legendre coefficients, and P n are the orthogonal Legendre polynomials. In this study, the Legendre coefficients, C n are derived by solving the linear system of equations at each equatorial pitch angle, α eq given by Equation 2. Note that the Legendre polynomials with even n-values are functions symmetric about an equatorial pitch angle of 90° and Legendre polynomials with odd n-values are antisymmetric functions (Chen et al., 2014). Since some of the highly asymmetric PADs likely result from errors in the flux measurements, all highly asymmetric PADs are excluded from this study. We use the approach outlined in Chen et al. (2014) to remove any highly asymmetric PAD fits. The asymmetry in the PAD fits is quantified by comparing each even C n coefficient with the consecutive odd coefficient C n+1 . If any even C n coefficient is less than the consecutive odd C n+1 , the resulting asymmetric PAD fit is excluded from this study. In Figure S2 in Supporting Information S1, an example of a highly symmetric PAD fit is shown, where each even C n coefficient has a greater magnitude than the consecutive odd coefficient. For comparison, Figure S2 in Supporting Information S1 also shows an example of an asymmetric PAD fit that is removed from the study, where some of the even C n coefficients have a smaller magnitude than the consecutive odd coefficients. Each of the PAD fits consists of Legendre polynomials and a coefficient up to order eight (n = 0, 1, 2, 3, 4, 5, 6, 7, and 8). Finally, we ensure that the PAD fits closely match the measured PAD values by removing any fits with fitted flux values over a factor of two higher or a factor of 10 lower than the measured flux values at the observed pitch angle. Examples of the PAD fits as well as the PAD data at the discrete equatorial pitch angles are illustrated in Figure 1. These PAD fits give a higher resolution characterization of how the PADs vary with the equatorial pitch angle than can be obtained by using the PAD data provided by either the MagEIS or REPT instruments themselves, which are only capable of measuring the flux at a limited number of local pitch angles.
In this paper the PADs fits are used to determine the values of the peak flux, and for the case of the butterfly PADs, the value of the local flux minimum. In addition, the fits are also used to determine the equatorial pitch angles where these flux peaks and local minimums occur. Note that since the odd Legendre polynomial coefficients, C n , used to produce the fits are not zero, the resulting PAD fits are not perfectly symmetric. Consequently, fitted butterfly PADs do not always have a local flux minimum at an equatorial pitch angle of precisely 90°. Similarly, the non-butterfly PADs may also have a peak flux that does not occur at an equatorial pitch angle of precisely 90°. However, due to the fits being composted of Legendre polynomials with even coefficients always greater than the odd coefficients, the resulting fits are highly symmetric with both the peak flux for nonbutterfly PADs and local flux minimum for butterfly PADs located at equatorial pitch angles within ±15° of 90°. 1. The flux data and the fits to the flux data both must have a local flux minimum at an equatorial pitch angle that is within ±15° of the 90° equatorial pitch angle. 2. The maximum flux in both the flux data and the fits to the flux data occur at a local maximum. 3. The local flux minimum in both the flux data and the fits to the flux data occurs at an equatorial pitch angle closer to 90° than the equatorial pitch angle of the flux maximum. 4. The flux data and the fits to the flux data both must have a normalized butterfly flux dip size, NDS > 0.05.
NDS is defined in Equation 3.
The normalized butterfly flux dip size, NDS, used to specify the depth of the flux dip that occurs in each of the butterfly PAD is defined as, Note that for butterfly PADs, where the local flux minimum occurs at a global flux minimum as illustrated in Figure 1d, NDS = 1. However, for butterfly PADs, where the local flux minimum does not occur at a global flux minimum as illustrated in Figure 1e, NDS<1.
This butterfly PAD selection criteria removes questionable butterfly PADs likely due to errors in the flux measurements (e.g., uncertainties associated with low count statistics). Figure 2 illustrates the number of non-butterfly PADs that satisfy the previous data selection criteria in Section 2.1 as well as the number of butterfly PADs that also satisfy the selection criteria outlined above in Section 3.1. Note that more butterfly and non-butterfly PADs are identified in the REPT data compared with that in the MagEIS data as illustrated in Figure 2. The level 3 REPT data provides the flux in 17 local pitch angle bins, which may contribute to the higher occurrence of PADs satisfying the data selection criteria in the REPT data compared to the level 3 data from MagEIS, which provides the flux data in only 11 local pitch angle bins. In addition, results presented in Baker et al. (2019 see Figure A3c) suggest that the flux below ∼10 3 /cm 2 /s/sr at 2-4 MeV energies is better resolved in the REPT data compared with the MagEIS data, due to the lower noise floor of the REPT instrument.

Two Discrete and Distinct Butterfly PADs
Most of the butterfly PADs identified in Figure 2 have a flux, which peaks at equatorial pitch angles of either ∼35° or ∼65°, similar to the PAD profiles shown in Figures 1c and 1d, respectively. Figure 3 shows the percentage of butterfly PADs at each energy, which have flux peaks at equatorial pitch angles ranging from 5° to 85° with a step size of 10° and a bin width of ±5°. In general, during geomagnetically quiet times on low L-shells, most of the butterfly PADs have a flux, which peaks at equatorial pitch angles of 65° ± 5°. Conversely, during more active times on higher L-shells, most of the butterfly PADs have a flux, which peaks at equatorial pitch angles of 35° ± 5° (see Figures 3m and 3n). Consequently, in the occurrence distributions of the butterfly PADs shown in Figure 3 there are two distinct peaks at equatorial pitch angles of ∼35° and ∼65°. These two distinct populations of butterfly PADs with flux peaks at ∼35° and ∼65° occur across all energies from 0.2 to 3.4 MeV. These results are consistent with the results presented in Ni et al. (2016), who used local pitch angle data from the REPT instrument, rather than the mapped equatorial pitch angle data used here, to show that at ≥1.8 MeV on high L-shells (L = 6), most butterfly PADs have a peak flux at local pitch angles of 37°-58° and at 58°-79° on lower L-shells (L = 4). The equatorial pitch angle where the flux reaches a peak value is identified in Figure 3 using the Legendre polynomials fits instead of using the local pitch angle data, since the PAD data from the REPT and MagEIS instruments is only available at a limited number of discrete pitch angles.

Spatial Distribution of Butterfly PADs
The butterfly PADs with a peak flux at equatorial pitch angles of ∼35° occur over a different L-shell and MLT range compared with those detected with a peak flux at equatorial pitch angles of ∼65°. Figures 4 and 5 illustrate the percentage of all butterfly and nonbutterfly PADs detected that are butterfly PADs with the peak flux at equatorial pitch angles of 35° ± 5° and 65° ± 5°, respectively. Here, the results are determined at 24 equally spaced MLT bins, with a bin with of 1 hr, and at 5 equally spaced L* bins, from L* = 3 to L* = 5, with a bin width of ±0.25. The percentage in a given MLT L-shell bin is defined as the number of butterfly PADs with a flux peak in the specified equatorial pitch angle bin, divided by the total number of butterfly and non-butterfly PADs detected in that MLT L-shell bin multiplied by 100. Only butterfly and nonbutterfly PADs that satisfy the data selection criteria outlined in Sections 2 and 3 are included. At most 100% of all PADs in a given MLT L-shell bin shown in Figures 4 and 5 are butterfly PADs with the peak flux at 35° ± 5° and 65° ± 5°, respectively.  Figure 4 illustrates that, in general, a higher percentage of the butterfly PADs with a peak flux at 35° ± 5° occurs at higher L* values and during more geomagnetically active times, consistent with the results presented in Figure 3. In addition, Figure 4 also shows that, in general, the highest percentage of these butterfly PADs are detected on the nightside, typically between 18 and 24 MLT. In general, a greater proportion of these butterfly PADs are also detected at higher energies. Figure 5 illustrates that the percentage of butterfly PADs detected with a peak flux at equatorial pitch angles of 65° ± 5° is much more evenly distributed in L* and MLT, compared to the results presented in Figure 4. A higher proportion of these butterfly PADs also occurs during more geomagnetically active times.

The Flux Dip Size of the Butterfly PADs
Butterfly PADs in the outer radiation belt can be created in three different ways. The flux at equatorial pitch angles near 90° can be eroded creating a dip in the flux as illustrated schematically by the top row of panels in  . Here, the total flux integrated over all equatorial pitch angle look directions gets smaller as the butterfly dip size increases as illustrated by the blue curve in Figure 6k.
Alternatively, butterfly PADs can also be created by processes, which preferentially enhance the flux of electrons with off-equatorial pitch angles as schematically illustrated in the middle row of panels in Figure 6. Here, this process results in a net increase in the total flux integrated over all equatorial pitch angle look directions as the butterfly PAD dip size increases as illustrated by the red curve in Figure 6k.
The schematic presented in Figure 6 indicates that by plotting the butterfly PAD dip size against the total flux integrated over all equatorial pitch angle look directions, it is possible to determine if in general larger butterfly dips are associated with lower total flux values, suggesting that the butterfly PADs are mostly created by loss due to magnetopause shadowing. Conversely, if the larger butterfly dips are generally associated with higher total flux values, then this indicates that the butterfly PADs are mostly created by a process that preferentially enhances the flux of the off-equatorial electrons, such as local acceleration by chorus and magnetosonic waves. In Figure 6, the butterfly PAD dip size, DS, and the flux integrated over all equatorial pitch angle look directions, J(α eq ) are derived using Equations 4 and 5.
In order to determine if the butterfly PADs are in general created by flux decreases as illustrated in Figure 6 (top row) or increases in flux as illustrated in Figure 6 (middle row), the dip size, DS, of each of the butterfly PADs is plotted as a function of the total flux integrated over all equatorial pitch angle look directions, J(α eq ).
Here, if flux enhancements are mostly due to processes that preferentially increase the flux of the electrons with off-equatorial pitch angles, then statistically most of the higher flux intensity butterfly PADs should be associated with higher butterfly dip sizes; see Figure 6 (middle and bottom rows). This type of butterfly PAD could be due to local acceleration of the off-equatorial electrons, if the acceleration process is faster than the pitch angle scattering processes, which reduce the butterfly PAD dip size, DS. Conversely, if lower total flux values are mostly associated with higher butterfly dip sizes, see Figure 6 (top and bottom rows), then these butterfly PADs cannot be due to processes that preferentially increase the flux of electrons with off-equatorial pitch angles. This type of butterfly PAD may be due to processes that erode the flux of electrons with equatorial pitch angles near 90°, such as drift-shell splitting combined with magnetopause shadowing.
The red dots in Figure 7 illustrate results for the butterfly PADs that have a peak flux at equatorial pitch angles of 35° ± 5°, and the green lines show least squares fits to these data points. Similarly, results for the butterfly PADs that have a peak flux at equatorial pitch angles of 65° ± 5° are illustrated by the blue dots, with the cyan line showing least squares fits to these data points. Note, the least squares fits are only determined when there are more than 10 data points so that the general trend of how the dip size, DS, varies with total flux, J(α eq ) can be detected. The results presented in Figure 7 illustrate that the butterfly PADs created with a flux peak at ∼65° (blue dots) have a much smaller dip size than those created with a flux peak at ∼35° (red dots). In addition, these results also show that both butterfly PADs have lower dip sizes generally occurring at higher total flux values, as indicated by the least squares fit trend lines, in Figure 7. Therefore, the results presented in Figure 7 suggest that the butterfly PADs are not generally due to a process that preferentially accelerates the off-equatorial electrons such as the acceleration process proposed in Xiao et al. (2015) (see Sections 5.1 and 5.2 for more details).

Correlation Between the Flux at Low and High Pitch Angles
The PADs presented in Figure 1 indicate that the flux can change sharply with the pitch angle as illustrated by the pancake PAD in Figure 1a, or remain almost constant with the pitch angle as illustrated by the flat-top PAD in Figure 1b. In Figure 8 and Figure 9, we examine how sharply the flux varies with the equatorial pitch angle for both the butterfly and non-butterfly PADs identified during the Van Allan Probe era. The gray squares in Figure 8 show the value of the flux of the non-butterfly PADs at equatorial pitch angles of 35° plotted against the value of the peak flux, at equatorial pitch angles near 90°. Similarly, the red squares show results for the butterfly PADs with a peak flux at 35° ± 5°, where the x-axes for these butterfly PADs indicate the flux value at the local flux minimum and the y-axes indicate the value of the peak flux. The solid line is simply y = x. By definition, the butterfly PADs all have points above the y = x line, since the flux at 35° ± 5° is always greater than the flux at the local flux minimum. Conversely, the non-butterfly PADs are all located below y = x, since these PADs all have a peak flux occurring near 90° that is greater than the flux at equatorial pitch angles of 35°. Points directly . Results are presented for L-shell bins of L* = 3, 4, and 5 ± 0.25 and energies of 0.2, 0.5, 1.1, 1.6, 2.1, 2.6, and 3.4 MeV. The red and blue dots represent butterfly PADs with maximum peak flux at equatorial pitch angles of 35° ± 5° and 65° ± 5°, respectively. The green and cyan lines illustrate least squares, LS, fits to the data, indicating that, in general, the greatest flux integrated over all look directions, J(α eq ), occurs for butterfly PADs that have the smallest flux dip size, DS. All of the results shown are for Dst = −30 ± 15 nT. Similar results for Dst = −60 ± 15 nT and Dst = 0 ± 15 nT as well as for L* = 3.5 ± 0.25 and 4.5 ± 0.25 are also presented in Figures S3, S4, S5, S6, and S7 in Supporting Information S1.  Figures S8, S9 S10, S11, and S12 in Supporting Information S1.  on y = x most likely indicate flat-top PADs with equal flux at equatorial pitch angles of 35° and 90°. The number of butterfly and non-butterfly PADs detected at each energy and in each L* bin is shown in red and gray, respectively. In general, Figure 8 illustrates that the flux at equatorial pitch angles of 35° and 90° are well correlated over all L-shells and energies as indicated by the trend of PADs with a high flux at 90° also having a high flux at 35°.  Consequently, panel (a) illustrates that at L* = 3 ± 0.25 and 0.2 MeV, a high percentage of PADs are of the pancake type where the flux at equatorial pitch angles of ∼90° is much higher than the flux at equatorial pitch angles of ∼35°. The red curves indicate similar results for the butterfly PADs shown in Figure 8. Here the x-axis indicates flux ratio bins, of the flux at equatorial pitch angles of ∼35° over the flux at equatorial pitch angles of ∼90°. Each panel also indicates the median value of these flux ratios for the butterfly (red) and nonbutterfly (blue) PADs, with the results shown in black representing the combined median flux ratio value of both the butterfly and non-butterfly PADs, identical to the median flux ratios shown in Figure 8. The flux ratios curves for the combined butterfly and nonbutterfly PADs are not shown and closely follow the nonbutterfly PAD ratio curves shown in blue. All results are for PADs obtained at Dst = −30 nT ± 15 nT. Similar, results for Dst = 0 nT ± 15 nT and Dst = − 60 nT ± 15 nT as well as for L* = 3.5 ± 0.25 and 4.5 ± 0.25 are also presented Figures S13, S14, S15, S16, and S17 in Supporting Information S1. The median flux ratios in Figure 8 indicate how close on average the flux is at equatorial pitch angles of 35° and 90°, for all the PADs detected in each individual L-shell bin and at each energy. These median flux ratio values are determined from the combined flux ratios from both the butterfly and non-butterfly PADs. For the butterfly PADs, each flux ratio is specified as the peak flux at 35° ± 5° divided by the value of the flux at the local flux minimum. For the non-butterfly PADs, each flux ratio is specified as the peak flux near 90° divided by the value of the flux at 35°. Interestingly, at L* = 3 ± 0.25 the flux at equatorial pitch angles of 35° and 90° is in general closer at higher energies than at lower energies. For example, Figure 8 shows that at L* = 3 ± 0.25, the 3.4 MeV electrons have a median flux ratio of 1.68, see Figure 8s, indicating that the 3.4 MeV electrons have on average a flatter PAD compared to the electrons at 0.2 MeV where the median flux ratio is 11.81; see Figure 8a. Conversely, on higher L-shells, it is lower energy PADs, which are flatter. For example, Figure 8 illustrates that at L* = 5 ± 0.25, the 0.2 MeV electrons have on average a flatter PAD with a median flux ratio of 1.56, see Figure 8c, compared with the higher energy 3.4 MeV electrons, which have a median flux ratio of 3.37; see Figure 8u.
The results presented in Figure 9 illustrate in more detail how the 35° over 90° equatorial pitch angle flux ratios change with energy and L-shell. The results in the left-hand column of Figure 9 show that at L* = 3 ± 0.25, the median flux ratios tend to decrease with increasing energy, with the largest drop occurring from 0.2 to 0.5 MeV.
These results indicate that at energies ≥0.5 MeV, the PADs are generally much flatter with similar flux values at equatorial pitch angles of 35° and 90°, consistent with the results presented in the left-hand column of Figure 8. Conversely, the right-hand column of Figure 9 illustrates that at L* = 5 ± 0.25, in general, both the butterfly and non-butterfly PADs gradually become flatter at lower energies as indicated by the continuous decrease in the median flux ratios with decreasing energy from 3.4 to 0.2 MeV. The PADs at L* = 4 ± 0.25, also generally tend to become slightly flatter at lower energies as indicated by the results in the middle column of Figure 9.

Discussion
The results presented in Sections 3 and 4 illustrate that the properties of butterfly and nonbutterfly equatorial PADs in the Earth's outer radiation belt vary with distinct properties as a function of L-shell, MLT, energy, and Dst. In this section, we discuss some of the possible physical processes that could be responsible for creating these electron PADs and further discuss possible reasons why their properties could vary with location, energy, and geomagnetic activity in the way observed.

The Butterfly PAD Population With a Flux Peak at 35°T
he results presented in Figure 3 show that most of the butterfly PADs detected during the Van Allen Probe era can be separated into two distinct populations with discrete occurrence peaks, those where the flux peaks at equatorial pitch angles near 35° and those where the flux peaks near 65°. These two distinct populations have different rates of occurrence as a function of L-shell and MLT, suggesting that each of these populations may be created by a different physical process. The butterfly PAD population with a peak flux near 35° has properties consistent with drift-shell splitting combined with magnetopause shadowing (West et al., 1973; see also Selesnick & Blake, 2002;Sibeck et al., 1987). Magnetopause shadowing and drift-shell splitting effects increase during more geomagnetically active times as the magnetopause moves inward, and the dayside magnetosphere becomes increasingly compressed. In addition, loss to the magnetopause will also be enhanced due to outward radial diffusion of the electrons to the magnetopause driven by ULF waves; this process also increases during more geomagnetically active times (Loto'aniu et al., 2010;Turner et al., 2012). Magnetopause shadowing, drift-shell splitting and outward radial diffusion should affect all energies of outer radiation belt particles. Consequently, the occurrence of this population of butterfly PADs at all energies, on higher L-shells, during more geomagnetically active times is consistent with these butterfly PADs being created by the combined effects of magnetopause shadowing, drift-shell splitting, and outward radial diffusion.
The results presented in Figure 4 also show that these butterfly PADs are predominantly detected on the nightside, generally between 18:00 MLT and 24:00 MLT. Note that since the asymmetric electron drift orbits pass closest to the Earth on the nightside, the Van Allen Probes spend more time on higher L-shells when their apogee is located on the nightside compared to when their apogee is located on the dayside. Consequently, the nightside occurrence peak in the butterfly PADs presented in Figure 4 is also consistent with butterfly PADs being created on higher L-shells. The results presented in Figures 3 and 4 are consistent with the results presented in Ni et al. (2016) also showing that more butterfly PADs are detected by the Van Allen Probes on the nightside at locations within 3° of the geomagnetic equator. These high L-shell nightside butterfly PADs presented in Ni et al. (2016) were identified from the data as having a peak flux somewhere between 37° and 58° at energies above 1.8 MeV. In this study, using fits to the PAD data, and a more rigorous selection criteria to better characterize the flux variation with pitch angle, we show that the majority of these butterfly PADs have a discrete and narrower peak flux at equatorial pitch angles in the range of 35° ± 5° at all energies from 0.2 to 3.4 MeV. The schematic shown in Figure 6 shows that loss processes, which erode the electron flux at equatorial pitch angles near 90°, such as drift-shell splitting in combination with magnetopause shadowing loss, will create butterfly PADs with greater dip sizes, DS, at lower values of the flux integrated over all equatorial pitch angles, J(α eq ). The trend of greater butterfly dip sizes, DS, at lower values of the flux integrated over all equatorial pitch angles, J(α eq ) is evident in Figure 7, consistent with drift-shell splitting in combination with magnetopause shadowing. This further supports our hypothesis that most of these butterfly PADs with a peak flux at 35° ± 5° are created by loss due to the combined effects of drift-shell splitting and magnetopause shadowing.

The Butterfly PAD Population With a Flux Peak at 65⁰
The occurrence of the other major population of butterfly PADs with a peak flux at 65° ± 5° is more evenly distributed across MLT and L-shells, with slightly more of these PADs being identified at higher energies as illustrated in Figure 4. Compared with the butterfly PADs with a peak flux at 35° ± 5°, a much higher proportion of these butterfly PADs occur on the dayside and on low L-shells even during low levels of geomagnetic activity. This suggests that they are not predominantly created by magnetopause shadowing loss. One possible explanation for this population of butterfly PADs with a peak flux at 65° ± 5° is that they are created by the combined effects of local acceleration by chorus waves as well as pitch angle scattering processes (see e.g., Albert et al., 2016). Local acceleration by chorus waves has been shown to be most efficient at equatorial pitch angles between 60° and 80° (see e.g., Albert et al., 2016;Horne et al., 2005;Thorne et al., 2013). Similarly, Horne et al. (2007) indicate that diffusive local acceleration of electrons from 0.1 to 3 MeV by magnetosonic waves would also be most efficient for electrons with equatorial pitch angles between 50° and 80°. Analysis of the PADs presented in Ni et al. (2020) further suggest that long lasting butterfly PADs are likely due to magnetosonic waves (see also Xiao et al., 2015). Consequently, these waves could be responsible for creating the population of butterfly PADs with a peak flux at 65° ± 5° (see, Horne et al., 2005). However, the observations presented in Figure 8 indicate that the dip size, DS, of these butterfly PADs tends to become greater at lower values of the flux integrated over all equatorial pitch angles, J(α eq ), inconsistent with an off-equatorial acceleration process that enhances the flux of 65° ± 5° relative to that at 90°; see Figures 6 and 7.
Here, we hypothesize that these electrons are locally accelerated by the chorus and magnetosonic waves, combined with pitch angle scattering processes due to hiss and chorus waves that act to reduce the size of the butterfly PAD dip, DS, so that at higher flux values, J(α eq ), the dip size, DS, of the butterfly PADs will be smaller. Several outer radiation belt studies, based on combined chorus wave momentum and pitch angle diffusion coefficients have also indicated that local acceleration by chorus waves in combination with pitch angle scattering processes may create PADs with small dip sizes and a peak flux near 65°, or even flat nonbutterfly PADs (see e.g., Albert et al., 2016;Horne et al., 2005;Jin et al., 2018;Li et al., 2014Li et al., , 2016Thorne et al., 2013).
While the majority of the butterfly PADs identified in this study with flux peaks at 35° ± 5° and 65° ± 5° have properties consistent with the processes discussed above, it is possible that some of these PADs are also created by other processes. For example, some of the butterfly PADs could be created by local acceleration processes that produce butterfly PADs with increasing butterfly dip sizes, DS, at higher integral flux values, J(α eq ) as illustrated in Figure 6 (middle and bottom rows). In addition, Artemyev et al. (2015) showed that butterfly PADs can be created by the current due to hot ions injected into the inner magnetosphere. This current of hot ions deforms the equatorial magnetic field and preferentially scatters electrons with the equatorial pitch angle near 90° down to lower pitch angles, as the particles' second adiabatic invariant is broken. The results presented in Artemyev et al. (2015) show that this process creates electron butterfly PADs on the nightside with energies from 1 to 3 MeV during geomagnetically quiet and active times. Such processes may result in the creation of butterfly PADs without causing significant overall loss or enhancement in the outer radiation belt's total flux integrated over all equatorial pitch angles, J(α eq ). Here, the butterfly dip size would be independent of J(α eq ). However, since the results presented in Figure 7 show that, in general, the butterfly dip size, DS, of the butterfly PADs tends to decrease with increasing total flux, J(α eq ), this suggests overall that the majority of the PADs are not created by either of these processes, and instead are likely associated with the combined effects of local acceleration and pitch angle scattering. Test particle simulations presented in Kamiya et al. (2018), show that a drift resonance interaction between long lasting monochromatic ULF waves and outer radiation belt electrons can in theory result in the formation of butterfly PADs that also could potentially contribute to some of the observed butterfly PAD presented in Figure 7.

The Correlation Between the Flux at Low and High Pitch Angles
The results presented in Section 4 illustrate that for both butterfly and nonbutterfly PADs the flux at low and high equatorial pitch angles remain highly correlated even as the flux varies by four orders of magnitude. This somewhat remarkable result indicates that the shape of the PADs remains very similar at flux levels that can vary by orders of magnitude. This result indicates that it may be possible to infer the flux of outer radiation belt electrons with equatorial pitch angles near 90° from measurements of the electron flux taken at lower pitch angles. One important implication and application of this result could be that outer radiation belt flux measurements taken by spacecraft in Low Earth Orbit (LEO) or Medium Earth Orbit (MEO), which are only able to detect electrons with low equatorial pitch angles, could potentially be used to estimate the flux of electrons at higher altitudes and at higher equatorial pitch angles. This approach has recently been tested in Claudepierre and O'Brien (2020); see also Pierrard et al. (2020Pierrard et al. ( , 2021. Similar results were also presented in Kanekal et al. (2001), who used 2 years of data collected from multiple spacecraft to show that the outer radiation belt flux observed by low inclination spacecraft at high altitudes in the equatorial plane is highly correlated with the flux on the same L-shell measured by high inclination spacecraft at much lower altitudes in LEO.

Why Does the Shape of the PADs Vary With Energy and L-Shell?
The results presented in Section 4 also indicate that the shapes of the butterfly and nonbutterfly PADs change clearly and characteristically with both energy and L-shells. Figures 8 and 9 show that in general, at L* = 3 ± 0.25, the flux values at low (∼35°) and high (∼90°) equatorial pitch angles are much closer at higher energies, indicating that on these low L-shells, the PADs are flatter between these pitch angles at higher energies. Conversely, Figures 8 and 9 also show that at L* = 5 ± 0.25, the flux values at low and high equatorial pitch angles are, in general, slightly further apart at higher energies, indicating that on these higher L-shells, it is the PADs at lower energies that are flatter. Note that when the butterfly PADs become flatter, the normalized butterfly dip size, NDS, is reduced. Here, we discuss the physical processes that may explain why the PADs flatten with increasing energy on low L-shells and with decreasing energy on high L-shells.
One possible explanation for the flattening of the PADs with increasing energy on low L-shells is that the electrons are scattered by plasmaspheric hiss waves with energy dependent pitch angle diffusion coefficients. Results presented in Malaspina et al. (2020) and in Ripoll et al. (2019Ripoll et al. ( , 2021 show that at low electron energies, ∼50-∼300 keV, the pitch angle diffusion coefficients due to plasmaspheric hiss waves are much higher at lower pitch angles than at higher pitch angles. Consequently, electrons on these low L-shells and with both low energies and low equatorial pitch angles will be preferentially scattered into the loss cone, eroding the flux of electrons with low pitch angles more than the flux of the equatorially mirroring electrons, creating PADs that have higher flux at equatorial pitch angles near 90°. However, at high electron energies, ≳1 MeV, the pitch angle scattering rates due to the plasmaspheric hiss waves are smaller and remain generally similar across all equatorial pitch angles (see e.g., Malaspina et al., 2020;Ripoll et al., 2019Ripoll et al., , 2021. Here, the electrons on the same low L-shells but with high energies, ≳1 MeV, will be weakly scattered by the plasmaspheric hiss waves across a wider range of pitch angles and are not able to create PADs with lower flux values at lower equatorial pitch angles even during long lasting intervals of plasmaspheric hiss scattering, so that these PADs remain flatter. The simulation results presented by Ripoll et al. (2019) support this explanation, showing that hiss waves in the plasmasphere can reproduce an observed top-hat PAD at 0.1 MeV, which has higher flux values at equatorial pitch angles near 90° compared with the flux at lower equatorial pitch angles as illustrated in Figure 1c. Conversely, at L*∼5, which generally corresponds to locations outside the plasmasphere during geomagnetic storms (see e.g., O'Brien & Moldwin, 2003), the results presented in Figures 8 and 9 suggest that the PADs tend to be slightly more peaked at higher energies. One possible explanation for these slightly more peaked PADs at higher energies is that on higher L-shells, the rate of pitch angle scattering is dominated by EMIC and chorus waves. Previous studies have shown that pitch angle diffusion rates due to chorus waves are higher for outer radiation belt electrons with lower energies (see e.g., Horne et al., 2005;Ripoll et al., 2021;Shprits et al., 2007), with similar high pitch angle diffusion rates across a wide range of equatorial pitch angles (see, Orlova & Shprits, 2010). At higher energies >1 MeV EMIC waves can play a dominant role in scattering the electrons into the loss cone (see e.g., Ukhorskiy et al., 2010). These EMIC waves have pitch angle diffusion coefficients much higher at lower equatorial pitch angles (see e.g., Ni et al., 2015;Ross et al., 2020), suggesting that these waves may contribute to the occurrence of top-hat PADs at energies >1 MeV, which have higher flux values at equatorial pitch angles near 90⁰. Consequently, on higher L-shells, typically outside of the plasmasphere region where the rate of pitch angle scattering is not affected by plasmaspheric hiss waves, the outer radiation belt PADs become increasingly more peaked at higher energies, especially at energies >1 MeV, consistent with the results presented in Figures 8 and 9.

Summary
In this paper, PAD data collected from the MagEIS and REPT instruments during the complete Van Allen Probes mission is used to characterize how the flux of outer radiation belt electrons depends on the particle's equatorial pitch angle, location, energy and geomagnetic activity. The results presented here mainly focus on the occurrence statistics and properties of butterfly PADs. In future work the properties of the non-butterfly PADs will be examined in more detail.
The PADs identified in this study have the following properties: 1. Butterfly PADs have a flux that has two distinct and discrete occurrence peaks at equatorial pitch angles of either ∼35° or ∼65°2 . The butterfly PAD populations with a peak flux near 35° mostly occur on the nightside at higher L-shells during geomagnetically active times with a dip size, DS, that is greater at lower total flux values, J(α eq ), consistent with the properties of butterfly PADs created by drift-shell splitting and magnetopause shadowing loss. 3. The butterfly PAD population with a peak flux at ∼65° occurs more uniformly across all L-shells and MLTs with properties suggesting these PADs are not predominantly created by magnetopause shadowing loss. These butterfly PADs may be caused by the combined effects of local acceleration, which increases the off-equatorial flux, and pitch angle scattering, which can scatter the increased flux to higher pitch angles, producing flux enhancements without continuously increasing the dip size of the butterfly PADs. 4. The flux at low and high pitch angles remains highly correlated even as the flux intensity varies by four orders of magnitude regardless of the shape of the PADs. 5. On low L-shells, for example, at L* = 3 ± 0.25, the PADs become increasingly more peaked at lower energies, consistent with the effects of the pitch angle and energy dependent pitch angle scattering rates resulting from plasmaspheric hiss waves. 6. On high L-shells, for example, L* = 5 ± 0.25, the PADs become increasingly more peaked at higher energies, especially at energies >1 MeV, consistent with loss dominated by EMIC waves, which preferentially scatter low pitch angle electrons at energies >1 MeV into the atmosphere.
In addition, our results also demonstrate that the shape of electron PADs remain statistically very similar even as the flux changes by up to four orders of magnitude, suggesting the utility of using energetic particle data from low earth orbit for characterization of the radiation belts at higher altitudes even in the heart of the Van Allen belts.

Data Availability Statement
All Supplementary Figures are available in the supporting material at the Zenodo data repository (https://doi.org/10.5281/zenodo.5899672).