Measurability of the Nonlinear Response of Electron Distribution Function to Chorus Emissions in the Earth's Radiation Belt
Abstract
We conduct test particle simulations to study the perturbations in a hot electron velocity distribution caused by a rising chorus element propagating parallel to the ambient magnetic field in the Earth's outer radiation belt. The wavefield is constructed from the nonlinear growth theory of chorus emissions of Omura (2021, https://doi.org/10.1186/s40623-021-01380-w), with additional considerations about saturation and propagation of the transverse resonant current being applied to model the subpacket structure. Using Liouville's theorem, we trace electrons back in time to reconstruct the evolution of electron velocity distribution at the magnetic equator. The electromagnetic hole created by nonlinear trapping and transport effects appears as a depression in the velocity distribution, aligned with the resonance velocity curve. We analyze the decrease of particle flux in this depression and estimate the energy resolution, pitch angle resolution, time resolution and geometric factor of particle analyzers needed to observe the perturbation. We conclude that particle detectors on current or recently operating spacecraft are always lacking in at least one of these parameters, which explains the missing direct observations of sharp phase space density depressions during chorus-electron nonlinear resonant interaction. However, with a dedicated experiment and appropriate measurement strategy, such observations are within the possibilities of the current technology. Similarity of the simulated density perturbation and a step function mathematical model is used to draw an analogy between the backward wave oscillator regime of chorus generation and the nonlinear growth theory.
Plain Language Summary
The plasma environment in the Earth's magnetosphere supports natural growth of various electromagnetic waves, including the whistler-mode chorus emissions, which consist of nonlinear chirping tones. These emissions can reach large amplitudes and play a major role in energization of radiation belt electrons. Nonlinear theories of chorus generation imply microscopic perturbations in resonant electron populations. A long-standing problem is that these predictions were never directly confirmed by experimental observations. Here, we analyze perturbations of electron distribution functions numerically, taking into account spacecraft measurements of short subpackets within each chirping element. We reveal distinct perturbations, which are just below the measurability limits of existing spacecraft instruments. We thus explain the current absence of direct measurements of nonlinear effects of chorus on the electron distribution functions. We also suggest measurement strategies for future spacecraft instruments that can increase the number of detected interaction events.
1 Introduction
The chorus electromagnetic emissions are a class of whistler-mode waves that occur naturally in planetary magnetospheres (Burtis & Helliwell, 1976; Hospodarsky et al., 2008; Tsurutani & Smith, 1974). They occur abundantly in the Earth's radiation belts, mostly in the midnight-to-noon sector (W. Li et al., 2009), and are recognized as an important driver of electron acceleration in this region (Baker et al., 2018; Thorne et al., 2013). They are responsible for the production of the relativistic “killer” electrons (Horne, 2007), as well as for electron precipitation into the upper atmosphere (Bortnik & Thorne, 2007; Ozaki et al., 2019), which creates weak, pulsating auroras (Ma et al., 2020; Miyoshi et al., 2015). The energetic electrons interact with counter-streaming whistler-mode waves through Doppler-shifted cyclotron resonance (Karpman & Shklyar, 1977; Tsurutani & Lakhina, 1997). The energy of the resonant particles ranges typically from about to for the lower frequency band (0.1–0.5 , where denotes the local electron gyrofrequency) (W. Li et al., 2010) and to for the upper frequency band (0.5–0.8 ) (Inan et al., 1992; Meredith et al., 2009). Between those two bands, a power gap is often observed in the chorus spectra (Burtis & Helliwell, 1969; Gao et al., 2019; Santolík et al., 2003), whose origin is possibly related to the formation of two anisotropic hot electron populations that generate each frequency band (Artemyev & Mourenas, 2020; J. Li et al., 2019; Ratcliffe & Watt, 2017).
The lower band part of the emissions often shows highly structured spectra consisting of coherent, high-amplitude discrete elements with rising or falling frequencies (Burtis & Helliwell, 1969; Santolík et al., 2004, 2009). Based on observational evidence and theoretical investigations, the chorus wave elements are known to be formed around the magnetic equator and to propagate bidirectionally to higher latitudes (Demekhov et al., 2020; Lauben et al., 2002; Santolík et al., 2004), approximately following the magnetic field lines (Hanzelka & Santolík, 2019; Santolík, Macúšová, et al., 2014). Observations of the Cluster spacecraft mission revealed a fine structure of chorus emissions with short wave subpackets within each element (Santolík et al., 2003). This fine structure has also been later observed by the Van Allen Probes (Crabtree et al., 2017; Santolík, Kletzing, et al., 2014).
Several theories have been proposed to explain the formation of chorus emissions. In the initial stage of the growth process, unstable anisotropic electron populations transfer energy to whistler-mode waves through the cyclotron resonance (Gołkowski et al., 2019; Kennel & Petschek, 1966; Omura et al., 2008). This mechanism is well described by the linear theory. In the second stage, whistler-mode waves with the highest amplitude start entrapping increasing numbers of resonant electrons, disturbing thus the gyrotropy (uniform distribution in the phase of the gyrating perpendicular velocity vector) of the electron distribution, which causes the formation of resonant currents (Nunn, 1974; Omura et al., 1991; Trakhtengerts et al., 2003). This nongyrotropic distribution carries a depression or elevation in phase space density, which is sometimes called an electromagnetic hole or hill (Omura et al., 2013). The resonant current modifies the dispersion properties of the wave and results in an enhanced, nonlinear amplitude growth, and in a drifting wave frequency. Simplified calculations of the resonant current have been carried out to obtain an estimate of the frequency drift and amplitude growth (Omura et al., 2008; Summers et al., 2012). Other approaches replace the electromagnetic hole with a step function in the parallel velocity distribution function, which then serves as the energy source for the nonlinear growth (Demekhov, 2011; Trakhtengerts, 1999). However, no theoretical approach so far has explained the formation of the subpackets within chorus elements. Observation-based models (Hanzelka et al., 2020) and kinetic particle simulations (Hikishima et al., 2010; Ke et al., 2017; Tao et al., 2017) can reproduce the fine structure, but the fundamental physics behind the subpackets is still unclear.
A mechanism of packet formation through wave beating has been proposed by Santolík (2008). This hypothesis is further supported by the analysis of Crabtree et al. (2017) who identified the simultaneous presence of multiple distinct waves within some chorus elements, especially toward the higher frequency end of the element. New statistical analyses of frequency variation and lengths of subpackets (Zhang, Mourenas, et al., 2020) further support the idea that the formation of fine structure can be at least partially explained by simple superposition of waves, especially in the case of short subpackets, with the longer ones being better explained by trapping-induced amplitude modulation (Tao et al., 2020). Nunn et al. (2021) conducted Vlasov hybrid simulations with two simultaneously occurring elements at similar frequencies (cf., Katoh & Omura, 2016; W. Li et al., 2011) and showed that the combined waveform features a rapid amplitude modulation, consistent with observations of the fine structure of chorus waves. It has been demonstrated with reduced kinetic theory (Mourenas et al., 2018) and by test particle simulations (Hiraga & Omura, 2020) that the presence of fine structure modifies the effectivity of electron acceleration, as the amplitude modulation and random phase jumps between subpackets disturb the trapped particle populations (Zhang, Agapitov, et al., 2020). Detailed understanding of the structure of chorus wave packets is therefore important for accurate predictions of energization and scattering of electrons in the Earth's radiation belts.
Particle analyzers on spacecraft orbiting in the Earth's magnetosphere have provided many new results in the recent years which improved our understanding of nonlinear wave-particle interactions. Fennell et al. (2014) have analyzed Van Allen Probes (RBSP) field and particle data and identified bursts of resonant electrons correlated with upper band chorus activity. Foster et al. (2017) combined theory and RBSP measurements of electron pitch angle distributions to estimate energization of relativistic electrons happening on time scales shorter than a typical chorus element. Perhaps the most direct observation of effects of nonlinear wave-particle interactions was provided by Kurita et al. (2018), who used 8-s measurements of electron fluxes provided by the Arase spacecraft during strong upper band chorus activity and observed transport of resonant electrons in the energy-pitch angle space, which followed the theoretical resonance curves. Nevertheless, spacecraft missions have been unable to detect depressions or elevation in the electron velocity distribution directly correlated with the nonlinear trapping within the duration of a single chorus element. In the related studies of nonlinear resonant interaction between protons and rising tone electromagnetic ion cyclotron (EMIC) waves, Shoji et al. (2017) found significant nongyrotropy in THEMIS proton data associated with occurrence of EMIC risers, which confirms the presence of an electromagnetic proton hole. However, to obtain sufficiently large number of particle counts, they had to average over time intervals containing several elements of the electromagnetic emission.
In this article, we study the perturbations to the hot electron distribution by means of test particle simulations under a fixed chorus wavefield. We include subpacket modulation of the wave amplitude in our model (Section 2) and track the time evolution of the distribution at the magnetic equator. Unlike most of the previous works on this subject, our focus is on obtaining full 2D energy-pitch angle distributions on subpacket timescales. In Sections 3.1 and 3.2, a stripe structure in the velocity distribution is observed and briefly described in agreement with previous kinetic simulations (e.g., Tao et al., 2017). In Section 3.3, we focus on the leading stripe, which is associated with the electromagnetic hole, and analyze the possibility to measure this perturbation with particle analyzers on recently operating spacecraft. Since the stripe in the phase space density has properties of a smoothed step function positioned near the resonance velocity, a successful observation of this feature would provide connection between the backward wave oscillator (BWO) theory of chorus growth, which relies on the existence of such a step function (Demekhov, 2017; Trakhtengerts, 1999), and the nonlinear growth theory of Omura (2021). We conclude that the pitch angle resolution of electron analyzers is insufficient to observe any of these stripes in their entirety. As it turns out, the most crucial deficiency comes from low particle counts on the short timescales of one subpacket, leading to the need of higher geometric factors of particle analyzers, or a change in measurement strategies.
2 Methods
2.1 Chorus Wave Field Model
To get meaningful results from test particle simulations, we need to construct a wavefield model that resembles chorus elements as observed by spacecraft. In Figure 1, we show an example of a rising lower band chorus element, which we used as a reference for our model. The spectrogram of the total magnetic field power spectral density in Figure 1a was obtained from the multicomponent continuous burst-mode measurements of the EMFISIS Waves instrument on the Van Allen Probes B spacecraft. In this mode, waveforms of all six electromagnetic field components were captured at a sampling rate of 35,000 samples per second. The waveforms were analyzed using 512-point Fast Fourier Transform with the Hann window, 50% overlap, and averaging over three neighboring spectra. Calibrated magnetic field waveforms were transformed into the coordinates linked to the background magnetic field (Figure 1b). They were subsequently band-pass filtered between 0.1 and 0.49 of the local electron cyclotron frequency. The instantaneous frequency (Figure 1c) was then obtained as the time derivative of the instantaneous phase of an analytic signal reconstructed using the Hilbert transform. Frequency data from the time intervals where the amplitude dropped below were cut off (Santolík, Kletzing, et al., 2014). The plasma frequency was obtained from measurements of the upper hybrid frequency (Kurth et al., 2015).
The example element is about long and rises from frequency of about to . After the initial long and weak subpacket, two well separated, high-amplitude packets appear, reaching about of the background field magnitude. With increasing wave frequency, the subpacket structure becomes less regular. The frequency follows a linear trend with irregularities appearing not only in between, but also inside subpackets, suggesting that this behavior cannot be caused purely by Hilbert transform applied on low-amplitude sections. The wave remains quasiparallel during the measurement, with wave normal angle at peak amplitudes starting at about and steadily decreasing.
The evolution of the frequency is a pure advection with group velocity , while the amplitude equation is modified by a nonlinear term on the right-hand side which includes the resonant current component parallel to the wave electric field. These equations are solved by the upwind method with a time step and a spatial step . On this grid, the Courant-Friedrichs-Lewy condition is satisfied over the whole domain while the numerical diffusion is kept reasonably slow (LeVeque, 2007).
The initial boundary values for each subpacket are given by the coupled set of modified chorus equations (Omura et al., 2009; Hanzelka et al., 2020). For detailed description of these equations, the notation and the complete set of input parameters, see Appendix A, where the results of Hanzelka et al. (2020) are briefly summarized.
To obtain the phase of the wave magnetic field vector, frequencies and wave numbers are integrated for each subpacket (we first integrate over time in the source, then over space away from the source). During particle simulations, these integrated phase angles are bilinearly interpolated at the position of each particle to obtain the angle between and . The phase difference between adjacent subpackets is not specified, as it changes with the movement of the source region.
The final model is presented in Figure 2 as plots of magnetic field amplitudes and wave frequencies. In the following particle simulations, we assume waves propagating in both directions from the magnetic equator, as shown in the amplitude plot in Figure 2b. Figures 2c and 2d show the frequency growth within the subpackets, ranging from to , where is the equatorial angular gyrofrequency. The simulation was stopped when the starting frequency of the next subpacket would become higher than ; this upper frequency limit corresponds to observations of lower band chorus with a spectral gap, but it is arbitrary in our simulation procedure as there is no natural cutoff for strictly parallel propagating whistler-mode waves. The simulated chorus element features a realistic subpacket structure with an irregular growth of the wave frequency (compare with Figure 1c) and with an upstream shift of the source region—these properties have been observed in numerical simulations (Hikishima et al., 2009) and also in spacecraft measurements (Demekhov et al., 2020; Foster et al., 2017; Santolík, Kletzing, et al., 2014).
2.2 Backward-in-Time Test Particle Simulation
To obtain a high-resolution velocity distribution at a given point in space and time, we use the backward-in-time simulation method (Nunn & Omura, 2015). In this approach, the equation of motion for electrons in an electromagnetic field is solved from a fixed final point (after the interaction with chorus) to an arbitrary initial point at which the wavefield disappears. The phase space densities corresponding to the velocity of electrons at the end of the simulation are mapped to the fixed final point with the use of Liouville's theorem (Goldstein et al., 2001): the velocity distribution function is constant along the trajectories of a Hamiltonian system. Therefore, if we know the initial, unperturbed distribution at all points along the magnetic field line and we choose the velocities of electrons so that we cover uniformly the velocity space at the final point, we can obtain the phase space density (PSD) of the electron distribution after the interaction with a chorus wave (see Figure 3 for a schematic explanation of the ideas behind the method). This method is obviously not self-consistent, as it requires a precalculated wavefield which does not respond to the changes in the electron velocity distribution. However, unlike in self-consistent kinetic simulations, we can choose to analyze an arbitrary section of the velocity space, thus increasing the resolution. And, maybe even more importantly, the initial distribution can be chosen a posteriori, meaning that we can obtain the perturbations for all possible distributions with just one simulation.
The test particles are advanced backward-in-time by a type of the Boris algorithm (Higuera & Cary, 2017) which does not use the common particle phase approximation (Zenitani & Umeda, 2018). The algorithm preserves phase space volume, which is an important property that allows us to safely apply Liouville's theorem in the numerical simulation. The effect of the mirror force on electrons is included through additional perpendicular components of the background magnetic field (Katoh & Omura, 2006). We also tested that a finer grid in the wave model has only negligible influence on the structure and magnitude of observed perturbations in the electron velocity distribution (peak value of relative PSD change in the simulation with one subpacket changed by less than ). However, some grid-dependent changes can be observed, especially at the edges of electromagnetic holes, as the chaotic trajectory of the resonant electrons strongly depends on their gyrophase at the onset of the resonant interaction.
The initial distribution in the particle phases is set to be uniform. stands here for the number density of hot electrons at the equator, (or ) and (or ) are the equatorial thermal velocities (or relativistic momenta) and is the magnetic dipole field strength along the field line. As noted by Kuzichev et al. (2019), this distribution has a different asymptotic behavior at relativistic velocities than an anisotropic Jüttner distribution. However, we are interested in particles at resonant velocities, which are only weakly relativistic. The loss cone is also not included in the model, because particle behavior at low pitch angles plays no significant role in the nonlinear growth theory.
The input parameters common to all simulation runs are as follows. The ratio of the plasma frequency to the equatorial gyrofrequency is (based on measurements shown in Figure 1) and the ratio of the density of hot electrons to the density of cold electrons is , which is a reasonable value for populations with thermal velocities of a few units to tens of keV (Juhász et al., 2019; Walsh et al., 2020). The parallel and perpendicular thermal relativistic velocities at the equator are and , resulting in an equatorial temperature anisotropy . The background magnetic field is assumed to be dipolar with equatorial strength at the Earth's surface of and the particles are propagating on the -shell with . The corresponding equatorial electron gyrofrequency evaluates to . These parameters are the same as in the wave model (see Appendix A).
All test particle simulations were done on a uniform grid in () space, where is the angle between the wave magnetic field vector and the perpendicular velocity of electrons. The parallel velocities were sampled in the range with 256 points, perpendicular velocities in the range with 256 points, and phases with 64 points (particles with initial velocities larger than are excluded from the sample). To confirm that this resolution is sufficient (Voitcu et al., 2012), we limited the range in velocities to focus on the resonance region of the first subpacket, and doubled the number of gyrophase grid points. Apart from some finer features associated with the phase mixing of lower and higher density regions, there was no noticeable difference. In Movies S1–S3, which depict the evolution of phase space hole, the parallel velocities are sampled on a much smaller range and with 1,024 points, the phases are sampled with 512 points and the perpendicular velocities are fixed to a single value. Each frame in these animations requires a new simulation run.
2.3 Particle Count Analysis
FPI | HOPE | PEACE HEEA | MEP-e | |
---|---|---|---|---|
MMS | RBSP | Cluster | Arase | |
Angular resolution (deg) | 11.25 | 18 | 15 (3.75) | 22.5 (3.5) |
Energy range (keV) | 0.01–30 | 0.001–50 | 0.01–26.4 | 7–87 |
Energy levels per sweep | 32 | 36 | 72 | 16 |
Time per energy level (ms) | 0.23 | 21 | 4 | 15.6 |
Geometric factor per bin () | 2.1 | 0.2 | 6.0 (1.5) | 0.7 |
- Note. The listed parameters are based on Pollock et al. (2016) for FPI (Funsten et al., 2013), for HOPE (Johnstone et al., 1997), for PEACE, and on (Kasahara et al., 2018) for medium-energy particle (MEP). The cylindrical MEP instrument has seemingly the best resolution, but it provides only 16 azimuthal channels with large dead zones, thus providing data from of each channel. For the rest of the instruments, the angular resolution refers to the polar angle bins.
The number of particles required to measure a significant decrease in phase space density can be defined in a variety of ways. Here we use the one-sigma (68%) confidence interval for a Poisson distribution approximated as , where is the number of particles measured on the perturbed distribution in a given energy-angle bin, and we require that the number of particles measured in the corresponding bin of the initial (unperturbed) distribution is at least by larger. Since relative change in number flux or PSD is directly proportional to the difference in number of measured particles, we get , where the subscript “0” denotes the initial state. The required number of particles measured in a single bin of the perturbed distribution can be then defined as . If we were to measure an increase in phase space density, we would relate the difference to the higher value of PSD, which would be the perturbed value, so we can generalize the requirement to .
We must also be mindful of the absolute number of measured particles. If we calculate the so-called exact confidence interval (Meeker et al., 2017), it will have different values for different , for example: for , for , and for . We can see that our approximated upper bound is noticeably underestimated for low , so we should always require in our assessment of measurability to make the approximations valid.
If we decided to work with a two-sigma (95%) confidence interval, the required number of particles would be about four times larger, with the upper bound of the exact confidence interval showing underestimation with nearly the same relative magnitude as in the one-sigma case.
3 Results
3.1 Response of the Electron Distribution to a Single Chorus Subpacket
We start by analyzing the evolution of the equatorial hot electron distribution during interaction with the first subpacket of our wave model from Section 2.1. To numerically simulate the wave-particle interactions, we use backward-in-time test particle simulations as described in Section 2.2. The final point of particle propagation (initial in the backward sense) is , , lying on the boundary between the first and the second subpacket. The time step is chosen to be a fixed value of . The particles remain near the equator where the magnetic field is almost constant, meaning that each electron gyroperiod is sampled with about 120 points.
During the resonant interaction, an electromagnetic hole forms in the space (see Movie S1 where the evolution of the electromagnetic hole in the first two subpackets is captured). Particles trapped in the hole are transported to lower parallel velocities, while resonant particles that flow around the hole increase their parallel velocity. For a highly anisotropic distribution, this translates to an increase of the pitch angle and the kinetic energy for trapped particles, and a decrease of and for untrapped particles (see Figures 4a and 4b for a sketch of the particle scattering and transport, and see Movies S2 and S3 for changes in the kinetic energy and pitch angle during the propagation of electrons through the first two subpackets). After integrating over particle phases, we obtain a 2D velocity distribution where the hole is not apparent anymore, and the perturbed part of the distribution appears as a stripe along the relativistic resonance velocity curve (see Figure 4c). The resonance velocity is defined as , where is the Lorentz factor and is the wave number. As expected from the inspection of the electromagnetic hole, the stripe consists of an increased phase space density part which is located close to resonance velocity corresponding to the initial frequency of the subpacket, and a decreased density part around resonance velocity corresponding to the wave frequency at the end of the subpacket (see Figure 4d). After integration over perpendicular velocities, we observe a smooth step-like feature in the reduced 1D distribution located around the resonance velocity (see Figure 4e).
3.2 Response of the Electron Distribution to a Full Chorus Element
In the analysis of the resonant interaction of electrons with a chorus element, we proceed in the same fashion as in the previous case of interaction with a single subpacket. Certain particles may now reach regions (or, in the forward time flow, come from regions) with significantly stronger magnetic field, so the time step must be appropriately smaller. A fixed value of was chosen, which translates to about 80 steps per electron gyroperiod at the boundary of the simulation domain where the magnetic field is about three times stronger than at the equator. Furthermore, some particles may reflect at their adiabatic mirror point and interact with the chorus element that propagates toward lower values.
We set the final point to , , which is right after the end of the last subpacket at the equator. The stripes created by the interaction with the subpackets are now distorted by adiabatic motion of particles that resonated with these subpackets farther from the equator (see Figure 5a; also see Movie S4 for subpacket-by-subpacket evolution of the perturbed distribution at the equator). In the 1D reduced distribution, the stripe-like perturbation related to single subpackets remains visible in those cases where the adiabatic motion of resonant electrons straightened the stripes and aligned them thus with the lines of constant parallel velocity (see Figure 5b). Subtraction of the initial distribution reveals an increase in phase space density at higher parallel velocities and a more prominent density decrease at lower velocities, especially around the resonance velocity corresponding to the last subpacket (see Figure 5d). Further analysis of the perturbations in the 2D distribution reveals that the stripes are overlapping, which makes them hard to distinguish, especially at resonance velocities corresponding to higher frequencies (see Figure 5c; also see the second half of Movie S1 to observe the mixing of phase space density between the first and second subpacket).
3.3 Analysis of the Required Resolution of Particle Measurements
The first step in this measurability analysis is to determine the time interval and the range of energies and pitch angles in which we observe large phase space density perturbations and, at the same time, large particle fluxes. If we assume the number of particles measured in each energy-angle bin to be Poisson distributed, then the statistical uncertainty of the measurement can be represented by the relative standard deviation . So, in order to make the measurement of a PSD perturbation significant at the one-sigma level, we need at least particles in each bin, where is the number flux (see Section 2 for a more detailed explanation of confidence intervals and the definition of number flux). This means that unlike in Figure 4d, we need to look at the relative perturbations of the phase space density and associated fluxes. The relative change in PSD after interaction with the first subpacket is shown in Figure 6a. The values of the number flux shown in Figure 6b, which are proportional to the particle counts, reach maximum at , . The quantity we use to determine the most suitable phase space region for measurements is . A plot of this quantity is presented in Figure 6c for the case of the first subpacket. By looking at the state of the perturbed distribution after interaction with each subpacket, we determined that the highest values of are found between subpackets 4 and 7. The case of the fifth subpacket is shown in Figure 6d. It needs to be stressed here that in our nonself-consistent simulation, we cannot determine how fast will the stripe on the left (associated with the first subpacket) decay, and so we should focus only on the newly formed generated right-hand stripe where the nonlinear interaction is taking place.
From Figure 6d we conclude that the optimal energy range for measurement of the phase space density decrease associated with growth of chorus subpackets starts at around and ends at about . In this energy range, the pitch angle width of the stripe is about 6°. A detector that always resolves the stripe should have a pitch angle resolution of half this width, which is 3°. Because the stripe does not lie on a constant pitch angle, the width of energy bins is also limited, but does not pose a major constraint.
In the second step of the measurability analysis, we compare the size and magnitude of the simulated perturbations with the resolution of an actual spacecraft instrument. Table 1 lists electrostatic particle detectors on several currently or recently operating spacecraft with science objectives related to wave-particle interactions in the magnetosphere. We can see that only the medium-energy particle-electrons (MEP-e) particle analyzer on the Arase spacecraft (Kasahara et al., 2018) has a suitable energy range, so we used this detector as a reference in further analysis. The MEP-e is a disklike electron spectrometer which filters incoming electrons by energy with an electrostatic analyzer and then uses avalanche photodiodes for particle detection. The time resolution of is close to the time duration of one chorus subpacket, meaning that we can use data from only one energy level to avoid smearing of the density stripe we wish to observe. With 16 logarithmically spaced energy levels from to , only the highest energy level is suitable to observe the structure depicted in Figure 6d (marked by a black arrow). The FPI detector on the MMS spacecraft and PEACE on the Cluster spacecraft provide better time resolution, but it compromises the particle counts they can obtain. For the purpose of comparison with our simulation, we assumed that the detector is oriented in such a way that it covers full in the pitch angle. Thus, the angular resolution in Table 1 is taken to be equivalent to pitch angle bins, although better resolution could be achieved with a different orientation and a more limited range of pitch angles. Changes of effectivity of the detector as a function of energy are not significant and are not considered here.
Since the measuring capabilities of electron analyzers on spacecraft are restricted to recording the particle count as a function of their kinetic energy and pitch angles, we will from now on inspect the perturbations of the electron distribution using these variables. Furthermore, all aforementioned simulation results represent the distribution at a single instant in time, while spacecraft particle analyzers typically require tens to hundreds of milliseconds to sweep through the relevant range of electron energies. To mimic the particle instrument behavior, we performed 16 backward particle simulations with the final time point being at the end of the fifth subpacket, and in each simulation, the final time point moved back by , adding up to a total time interval of . The chosen energy bin from to was covered with 16 logarithmically spaced grid points and we used 256 linearly spaced points for the pitch angle quadrant (64 points per single MEP-e angular sector) and 128 linearly spaced points for particle velocity phases.
Unfortunately, due to the large dead zones of the azimuthal channels of the MEP-e analyzer, it is immediately clear that we cannot resolve the perturbation in flux, because we have only one angular bin per a channel (see Table 1). Therefore, we had to base our simulation on a hypothetical detector which preserves the parameters of MEP-e, but has more azimuthal channels and negligible dead zones. We proceeded to refine the pitch angle bins to find out a compromise between high particle counts and good angular resolution. The ratio fluctuated as we were refining the angular grid, depending on whether the energy-angle bin fell into the center of the stripe or on the edge. However, with bin size or smaller, the ratio was growing steadily, signifying that by further resizing the bins, we are only losing particle counts and not improving the resolution of the stripe anymore. For bins of size between and , the ratio fluctuated between 1 and 2, with expected number of particles falling down to 3 or 4 per bin in the region of the PSD decrease. To achieve reasonably large particle counts, and also to prevent errors in estimation of confidence intervals, we further had to increase the geometric factor 10 times with respect to the original pitch angle resolution (i.e., the effective aperture of the detector has been increased 10 times), making it comparable to the PEACE detector.
In Figures 7a–7c, we can see an example with angular bin size of . A decrease in phase space density is observed between , followed by fluctuations at larger values of pitch angle. Thus, with the increased geometric factor and angular resolution, the phase space density stripe can be resolved, and the particle counts shown in Figure 7b are above 50. As shown in Figure 7c, the required particle counts for one-sigma significance are met at the lowest point of the dip. We can reach even , corresponding to a two-sigma significance. The points where the nonlinear wave growth is taking place are marked by green diamonds in the plot.
To get a more complete picture of the stripe as it could be seen by particle detectors, we further assumed four detectors operating simultaneously at four logarithmically spaced energy levels from to . We preserved the accumulation time step and energy level spacing as given in the MEP-e specifications, and we used the improved pitch angle resolution of and the increased geometric factor of as in Figures 7a-7c above. As expected from the above analysis, the PSD decrease is well resolved in both energies and pitch angles and the ratio reaches values much smaller than one at each energy level.
We conclude that with the particle detectors mounted on currently or recently operating spacecraft represented by Arase, Van Allen Probes, Cluster and MMS, a direct measurement of perturbations in electron velocity distribution due to interaction with a chorus subpackets is not possible, mainly due to insufficient angular resolution and low upper energy limits. With a dedicated experiment, however, an in situ observation of these perturbations should be within the current technical possibilities. Such an experiment would consist of particle analyzers measuring electrons with energies up to at least with a polar angle resolution better than about , and with a geometric factor of at least (i.e., 10-times better than MEP-e). Unlike with FPI on MMS, obtaining a fast 3D distribution would not be of essence: during high chorus wave activity, the analyzer would operate in a fixed high energy mode with accumulation times comparable to the theoretical time duration of a single subpacket. We believe that with this type of experiment, it should be possible to obtain direct confirmation of the processes described by the nonlinear growth theory of chorus rising tone emissions.
4 Summary and Discussion
We have studied the nonlinear interactions between a rising-frequency chorus element with a subpacket structure and a hot electron population. Each subpacket traps resonant electrons from a low phase space density part of the anisotropic velocity distribution and transports them to higher energies and pitch angles, while untrapped resonant particles experience scattering to lower energies and pitch angles. In result, we observe stripes of increased and decreased phase space density forming along the elliptical resonance curve, with increased PSD at the resonance corresponding to the initial frequency, and decreased PSD at the resonance corresponding to the terminal frequency of the subpacket. With each subsequent subpacket, some of the electrons from the phase space density depletion become trapped again. Due to this consecutive trapping, a very prominent depletion is formed at the position of the leading subpacket as the system evolves. Since we are observing the evolution from the equator, the perturbations which formed at higher latitudes will arrive to the observation point with adiabatic correction to pitch angle. Therefore, the stripes associated with previous subpackets do not follow resonance curves anymore, with their high pitch angle parts being shifted toward larger parallel velocities.
We further considered the possibility to measure the prominent PSD depletion by spacecraft particle instruments. Such measurements require energy and pitch angle range which covers the part of the depletion where the maximum change to particle flux is expected. Moreover, accumulation times are limited to a window in the order of tens of milliseconds due to the fast propagation of the perturbation through the velocity space. It turns out that none of the available instruments meets the criteria for a statistically significant observation. This negative result can be seen as an explanation for the lack of direct observational evidence of these perturbations. Nevertheless, we showed that with the energy resolution of near the energy level of (similar to the MEP-e instrument on Arase), a pitch angle resolution of (similar to MEP-e, but without dead zones), and with a geometric factor of (similar to that of the PEACE HEEA instrument on Cluster), it would be possible to achieve a conclusive measurement results on the phase space density depletion.
In our test particle simulation, we have used a fixed model of the wavefield of the chorus element. Such nonself-consistent approach can yield a realistic evolution of the electron velocity distribution function only if the model well represents the actual chorus wave fields present in the radiation belts. As it is technically impossible to directly measure a large time-space section of the wavefield, we can make a comparison only with kinetic simulations. Looking at some recent particle-in-cell (PIC) simulations (Hikishima et al., 2010; Tao et al., 2017), we can confirm that the phase space density decrease related to the last subpacket is well apparent in their results, while the rest of the stripe structure is only marginally visible. Our hypothesis here is that the stripes might decay while generating a broadband whistler-mode spectrum through a wave growth mechanism driven by increased pitch angle anisotropies along the stripes. Because this part of the spectrum is not present in our wave model, the additional stripes cannot decay in the test particle simulation. Based on the proposed hypothesis and the available PIC simulation results, we believe that observable effects of chorus subpackets on electron distributions primarily occur on the leading stripe with a large phase space density decrease. We note here that due to the sequential nature of our wave model and the drifting source assumption, each subpacket represents a rising-frequency plane wave with its phase independent of the previous packet, and so the limiting effect of random phase jumps on nonlinear electron acceleration (Zhang, Agapitov, et al., 2020) is naturally included. It should be also mentioned that the multi-wave structure (Crabtree et al., 2017) of chorus elements leads to higher variability in amplitude and packet length than assumed in our simple model, which could result in even more pronounced overlaps between resonance widths of each subpacket and add further smearing of the stripe structure. Unfortunately, the lack of observational evidence for the stripe structure cannot be used as evidence for multi-wave models, since this structure shows even smaller perturbation amplitudes than the PSD decrease analyzed in Section 3.3.
The presented analysis of measurability included some simplifications, mainly by assuming a special orientation of the electron analyzers with respect to the ambient magnetic field. A different inclination of the spectrometer could limit the total range of pitch angles by removing intervals around and . This would be acceptable since the interesting part of the distribution lies approximately within from the perpendicular direction. Also, we referenced mainly electrostatic analyzers, but given the required energy range and the advantage of a fixed energy operation mode, semiconductor detectors similar to RAPID (Wilken et al., 2001) on Cluster or FEEPS (Blake et al., 2016) on MMS might be more suitable for this task. However, the software on both RAPID and FEEPS was designed to provide data in time steps which are orders of magnitude above the subpacket time scale, so they were not considered in this paper.
We further note that the wavefield we used was based on a chorus event observed in region with a very low plasma-to-gyrofrequency ratio , and a very large perpendicular thermal velocity was assumed as well. These choices increase the resonance velocity and shift fluxes toward higher energies. A different wavefield with and was tested, which resulted in the optimal values starting at about , which is still outside the range of most of the available electrostatic analyzers.
Another simplification we used is the limitation to a single model of hot plasma distribution. Based on previous statistical results (W. Li et al., 2010), the phase space density of 10 – electrons measured between and near the equator ranges from about to during geomagnetically disturbed conditions. In our model, the average density in this energy range is about at , and so assuming even larger densities to achieve higher particle counts would make our estimates unrealistic. Further concerning the hot plasma distribution, we assume that it is smooth in its initial state. Given that there is a continuous presence of whistler waves in the radiation belts, occurrence of some fluctuations in the phase space density should be expected. The coherent nature of chorus will always lead to organization of the electrons and to creation of the electromagnetic hole, giving us confidence that the leading stripe structure will still be prominent, but the stripes outside of the active growth region might be further distorted by the ambient fluctuations.
The obtained phase space density decrease at the leading subpacket is similar to the step-like perturbation assumed in the BWO (Backward Wave Oscillator) chorus growth theory (Demekhov, 2017). However, the step in the 1D distribution (see the inset plot in Figure 5b) is not very sharp due to the relativistic effects included in the resonance velocity. Relativistic versions of the BWO theory formulas would be needed to make a comparison with our results. As this phase space density perturbation should be observable with a dedicated particle experiment, the same experiment should be also able to confirm experimentally the existence of the (relativistic) step function from the BWO theory. Positive outcome from the experiment would show a clear relation between the step function paradigm and the electromagnetic hole approach. The recent work of Zonca et al. (2021) presents a general theoretical framework for chorus excitation which hints at a fundamental connection between gyrotron backwave oscillators and other interpretations of chorus chirping (Omura & Nunn, 2011; Vomvoridis et al., 1982), supporting thus further experimental and theoretical research in this area.
In conclusion, our simulations show that the electron distribution function carries an imprint of the fine structure of nonlinear wave packets within chirping chorus elements. Spacecraft observation of the full stripe structure appearing in our simulated distribution would complement the research of subpacket formation in chorus elements by providing information about resonant currents, but unfortunately, such measurements are not feasible. On the other hand, the major perturbation associated with the front of the element provides connection between the BWO theory and the nonlinear growth theory of Omura (2021), and the possible detection of such perturbation is just beyond the possibilities of the currently available instrumentation. Better geometric factors of energetic electron analyzers are needed, along with operational modes focused on gathering large particle counts in relevant energy and pitch angle ranges. These improvements are thus needed to further support the validity of current theories of the chorus growth mechanism.
Notation
-
- equatorial temperature anisotropy.
-
- amplitude of wave magnetic field.
-
- speed of light in vacuum.
-
- elementary charge.
-
- position of the source of the -th subpacket along the field line.
-
- resonant current component parallel with the wave electric field.
-
- resonant current component parallel with the wave magnetic field.
-
- computed for .
-
- computed for .
-
- wave number.
-
- electron mass.
-
- hot electron number density.
-
- depth of the electron phase space hole.
-
- factors entering the calculation of the inhomogeneity .
-
- inhomogeneity ratio (Omura et al., 2008).
-
- negative value of the inhomogeneity ratio in the source.
-
- parallel relativistic particle momentum .
-
- perpendicular relativistic particle momentum .
-
- parallel thermal relativistic velocity, with subscript “eq” denoting the equatorial value.
-
- perpendicular thermal relativistic velocity, with subscript “eq” denoting the equatorial value.
-
- group velocity of whistler-mode wave.
-
- cyclotron resonance velocity.
-
- typical perpendicular velocity of particles, appears in the nonlinear growth theory in the function which replaces the perpendicular factor of the electron velocity distribution function.
-
- Lorentz factor.
-
- Lorentz factor for .
-
- nonlinear growth rate.
-
- difference between wave magnetic field phase and perpendicular particle velocity phase.
-
- vacuum permeability.
-
- wave angular frequency.
-
- normalized wave amplitude .
Acknowledgments
This work has received funding from the European Union's Horizon 2020 research and innovation programme under grant agreement No. 870452 (PAGER). M. Hanzelka, I. Kolmasova, and O. Santolik further acknowledge support from the Czech Academy of Sciences through the Mobility Plus grant JSPS-19-05, and from the MEYS grant LTAUSA17070. M. Hanzelka and O. Santolik acknowledge support from Charles University through the GA UK project No. 64120. Y. Omura acknowledges support from the JSPS KAKENHI Grant 17H06140.
Appendix A: Chorus Wave Model Equations
The average momentum at distance is related to the thermal momentum as . The distribution in Equation A9 is the same one as in our particle simulations.
Finally, the evolution equations are solved by the upwind method with time step and spatial step . The equatorial strength of the dipole field at the surface of the Earth is chosen as , L-value of the field line is . Further parameters are chosen as follows: , , , , , , .
Open Research
Data Availability Statement
The Van Allen Probe data are publicly available from the NASA's Space Physics Data Facility, repository https://spdf.gsfc.nasa.gov/pub/data/rbsp/. Data obtained from numerical simulations are available for download at https://doi.org/10.6084/m9.figshare.14315936.v1.