1 Introduction

Chorus emissions are coherent electromagnetic waves propagating in the whistler mode which are frequently observed in the inner magnetosphere, typically in the range of L shells from 4 to 8 (Kasahara et al., 2009; Santolík et al., 2003; Tsurutani & Smith, 1974). They can induce both acceleration and losses of energetic electrons in the radiation belts (Tsurutani et al., 2009; Turner et al., 2013) through nonlinear interactions (Summers et al., 2013). These processes are sensitive to the frequency-time structure of the chorus wave packets (Tao et al., 2013), which therefore needs to be well understood in order to fully comprehend the dynamics of the radiation belts. The fine structure of chorus elements which rise in frequency has been discovered from high-resolution measurements of the Cluster spacecraft (Santolík et al., 2003, 2004) which show that each element of the discrete emission consists of several subpackets with growing wave frequencies. The subpacket structure of chorus has been confirmed by recent analyses of multicomponent measurements of chorus by Van Allen Probes (Foster et al., 2017; Omura et al., 2019; Santolík et al., 2014). This fine structure has also been observed in full particle simulations (Hikishima et al., 2009, 2010, 2017) and hybrid simulations (Katoh & Omura, 2016). A feature unique to the simulations, not yet observed by any spacecraft missions, is the movement of the source to the region upstream of the wave, which happens along the frequency growth.

To explain the features of chorus emissions discovered in numerical simulations and spacecraft measurements, the nonlinear growth theory has been developed (Omura et al., 2008, 2009). This theory recognizes the inhomogeneity of magnetic field along a field line as the main controlling factor for the formation of an electromagnetic electron hole in the velocity phase space, near the resonance velocity (for evidence of the formation of the hole, see the numerical simulations of, e.g., Harid et al., 2014; Hikishima & Omura, 2012). Resonant electrons, gyrophase-bunched in the electron hole, produce a resonant current which causes the amplitude and frequency growth of the whistler mode wave. The nonlinear growth theory gives values of frequency sweep rates and amplitudes of chorus elements which are in good agreement with in situ observations (Foster et al., 2017; Kurita et al., 2012; Yagitani et al., 2014). It has also been applied to explain the fine structure of electromagnetic ion cyclotron (EMIC) emissions, which, similarly to chorus, consist of several subpackets (Nakamura et al., 2015; Omura et al., 2010). The subpacket structure of EMIC waves was analyzed numerically by Shoji and Omura (2013) and they also presented an idea that the subpackets could be produced by a repeated triggering process caused by the radiation from phase-organized protons which are continuously being released from the interaction region.

In the present study we use the nonlinear growth theory to develop a simple model of the fine structure of rising tone chorus emission, taking inspiration from the idea of subpacket formation in EMIC waves presented by Shoji and Omura (2013). The evolution of the wave amplitude and wave frequency inside a single subpacket in the source region is described by the so-called chorus equations, derived by Omura et al. (2009). Wave propagation and convective growth is modeled with advection equations. The fundamental assumption employed in the present model is that the resonant current, produced through wave-particle interaction, carries the information about the wave vector and frequency of the emission and can act as a helical antenna and radiate a new coherent wave during their upstream propagation. The resonant current appears as a nonlinear term on the right-hand side of the advection equation for wave amplitude. Similar idea (i.e., the resonant current acting as an antenna) already appeared in the seminal paper of Helliwell (1967), but they did not connect it with the nonlinear growth theory, which was not yet fully developed at that time. Trakhtengerts et al. (2003) analyzed the frequency shift due to this antenna effect and estimated the amplitude of the emitted radiation; however, they did not consider it as a possible cause for the subpacket structure. Here, some further assumptions are made to separate the newly radiated wave from the previous subpacket, and the optimum amplitude derived by Omura and Nunn (2011) is used to introduce saturation effects into the model. Chorus elements obtained from the numerical solution show that between adjacent subpacket, there are small, local drops in the otherwise growing frequency, which is a feature that seems to be also indicated by the measurements of the Van Allen Probes (Foster et al., 2017; Santolík et al., 2014). The upstream shift of the source region, previously obtained in some full-particle simulations (Hikishima et al., 2010; Ke et al., 2017), is also present in the model.

This new model of the subpacket structure of the chorus emission is introduced in section 2, which is further divided into three subsections that deal with the evolution equations for chorus, the resonant current, and the proposed sequence of processes that occur during the growth of a chorus element. In section 3 we present the numerical solution of the differential equations describing the new model, focusing on its unique features, namely the movement of the source region to the upstream and the inversion of frequency growth between subpackets. Section 4 is dedicated to the comparison of the modeled chorus element with Van Allen Probes observations of rising tone chorus emissions in the radiation belts. In section 5 we further discuss the advantages and shortcomings of the presented model and conclude our main results.

2 Model of a Chorus Element

2.1 The Evolution Equations

We are studying the evolution of wave frequency ω(h,t) and wave amplitude B_w(h,t) of a coherent electromagnetic whistler mode wave propagating parallel to a background dipole magnetic field through a one-component plasma with a constant number density of electrons. Distance h is measured along a magnetic field line, starting at the equator, t is the time. Following Summers et al. (2012), we describe the evolution with two coupled advection equations:

(1)

(2)

where V_g is the group velocity of a whistler mode wave, μ₀ is the permeability of vacuum, and J_E is the resonant current density component parallel to the wave electric field. The first equation simply states that the frequency is constant in a frame of reference moving with the group velocity, which is a consequence of the ray approximation (Lighthill, 1965). A detailed derivation of the second equation has been given by, for example, Nunn (1974) or Omura et al. (2008). Following Foster et al. (2017), we use Equations 1 and 2 to describe the evolution of a single subpacket, not the whole chorus element, which was done in previous studies, for example, Summers et al. (2012).

The time evolution of B_w and ω in the source is given by the chorus equations of Omura et al. (2009). To obtain the equation for ω, we start from the definition of the inhomogeneity ratio

(3)

where Ω_w is the normalized wave amplitude defined by Ω_w = eB_w/m_e, e denotes the elementary charge, m_e denotes the electron rest mass, and c is the speed of light in vacuum. The explicit forms of parameters s₀, s₁, and s₂ are given in Omura et al. (2009), Equations 11–13-11–13. Further, we will assume a parabolic approximation of the magnetic field strength along field lines, allowing us to define the dependence of electron gyrofrequency on the distance along field line as

(4)

where Ω_e0 = eB_eq/m_e is the equatorial electron gyrofrequency, B_eq is the magnetic field strength at the equator, and a comes from the small-latitude Taylor expansion of the magnetic field and is given by a = 4.5/(LR_E)², with R_E being the Earth's radius. Consequently,

(5)

We will require that |J_E| is maximized in the source, which is located in the distance h_i, where i indexes the subpackets. The maximum of |J_E| is achieved with (Omura et al., 2008) S ≈ −0.41 ≡ −S_max. We can now substitute this value on the left-hand side of Equation 3 to obtain, using also Equation 5, the first chorus equation

(6)

The second term on the right-hand side is not present in the derivation of similar equation presented in Omura et al. (2009), because in Equation 6 we have allowed the source to be located away from the equator.

The second chorus equation uses the concept of the threshold amplitude, which remains unchanged for h_i ≠ 0, so we can write (Omura et al., 2009)

(7)

Here Γ_N represents the growth rate defined by

(8)

As we will show in the next subsection, Γ_N depends on both Ω_w and ω, which causes a strongly nonlinear growth.

2.2 Resonant Current

The interaction between resonant electrons and whistler mode waves leads to the depletion of trapped electrons from the phase space, which is sometimes called the electromagnetic electron hole (Omura & Summers, 2006). The hole introduces a nongyrotropy into the velocity distribution, which then manifests through the appearance of the resonant current density J_R. It is useful to decompose this current density into the components J_E and J_B which are parallel to the wave electric and magnetic fields, respectively. The J_E component is connected to the growth of wave amplitude, as we have seen in Equation 2, and J_B causes the growth of wave frequency. They may be expressed as (Omura et al., 2008)

(9)

(10)

where ζ is the gyrophase angle defined with respect to the wave magnetic field, and ζ₁(S), ζ₂(S) set the left and right boundaries of the separatrix in the v_‖(ζ) phase portrait. The quantity J₀ depends on the distribution of hot electrons trapped by the wave. Here we follow Summers et al. (2012) and assume a fully adiabatic evolution of a hot electron distribution, chosen to be bi-Maxwellian in momenta, to define

(11)

where

(12)

carries information about the distribution function and

(13)

is the equatorial anisotropy of the hot electron distribution. The other quantities we introduced in Equations 11 and 12 are as follows: average perpendicular electron velocity V_⊥0, wave number k, resonance velocity V_R, Lorentz factor γ_R of an electron propagating with the resonance velocity, dimensionless parameter χ² = 1 − ω²/c²k² = 1 − 1/n² (where n is the refractive index of a whistler mode wave), number density N_he of the hot electron population, depth of the electron hole Q, equatorial perpendicular thermal velocity U_{th, ⊥ eq}, and equatorial parallel thermal velocity U_th, ‖eq. The wave number of a parallel whistler mode wave in cold plasma can be approximated as (Stix, 1992)

(14)

As a consequence of Equations 7, 9, 11-7, 9, 11, and 14, the nonlinear growth rate Γ_N defined in Equation 8 can be written explicitly as

(15)

showing a direct proportionality to . The constant J_E,max ≈ 0.98 gives the value of J_E = −J_E,maxJ₀ at S = −S_max and can be obtained by numerically evaluating Equation 9.

The particles which interact with the whistler wave have velocities and gyrophases that match the first-order resonance condition for a wave whose spatiotemporal structure is given by ω and k. Therefore, the particle bunches (and the depletion created by the bunching) form a helical shape in space on which are imprinted the wave frequency and wave vector of the interacting wave. Such helix acts as an antenna radiating a right-hand circularly polarized wave on this frequency. The use of helical antennas for creation of circularly polarized electromagnetic signals is a well-known concept in radio science, proposed in the 1940s by Kraus (1949). To get an estimate on the strength of the electromagnetic field radiated from the antenna, we will follow Yagitani et al. (1992) who computed the electric field of L-mode and R-mode plasma waves radiated from a current sheet on the background of a homogeneous magnetic field. Focusing on the R-mode, we can rewrite the result of Yagitani et al. (1992) as

(16)

This is the response of the electric field to a current distribution given by , where δ(z) is a delta distribution with units of inverse length and has the units of current density times length. Since we are not interested in the direction of the electric field vector, we have simplified the formula by assuming that J_s points along the x axis, leading to E_δ(z) having only one nonzero component in our coordinate system. To obtain the field radiated by the helical resonant current, we just have to realize that the electric field, E(z), will always point in the direction of the current at each point along the z axis, which coincides with the helical axis. Therefore, we only need to substitute the δ distribution with a more realistic distribution of the magnitude of the current. With the resonant current distribution given as

(17)

(18)

with σ_J being a characteristic width of the distribution, we can obtain the total radiated field at a point z → ∞ (far enough from the antenna) by integrating over the current distribution,

(19)

And since we have formulated the evolution equations in the terms of wave magnetic field, we can now use the relation c|B_tot|/n = |E_tot| to obtain

(20)

The quantity J_peak represents the peak value of the current density distribution, which may be obtained from a numerical simulation.

With a uniform distribution of the current

(21)

we would get

(22)

The strength of the magnetic field of the emitted wave is directly proportional to the length of the helix. This is in agreement with the strength of electromagnetic field of circularly polarized waves radiated from a helical antenna as derived by Kraus (1949), Equation 27.

We should note that the nonlinear advection Equation 2 already describes the amplitude of the whistler wave amplified by a (helical) resonant current, that is, it captures the antenna effect. As the particles which constitute the current propagate upstream of the wave, they should radiate continuously, spreading thus the source region of the chorus emission. However, the theory described in section 2.1 requires a point source. Therefore, we have to discretize the motion of the source and treat the current independently when it is released from the sharp amplitude gradient at h_i. This is discussed in the next section.

2.3 Model of the Subpacket Structure

We envision the formation of the subpacket structure of the whistler mode chorus as follows. Initially, the electromagnetic emissions in the equatorial region are dominated by incoherent noise. Through interaction with hot electrons, the amplitude of the noise grows according to the linear growth theory with a rate γ_L, which maximizes at the equator, as it was shown by numerical simulations (Hikishima et al., 2009; Katoh & Omura, 2016). After a certain time the linear growth produces a coherent emission with a wave amplitude that reaches the threshold amplitude (Omura et al., 2009)

(23)

where ω_phe denotes the plasma frequency of hot electrons. Ω_w > Ω_thr expresses the necessary condition to start the nonlinear growth rate stage—below this threshold value, Equations 6 and 7are not valid. Initially, ∂ω/∂t = 0 and ∂Ω_e/∂h = 0 at the equator, then S = 0 as a consequence of Equation 3. Under such conditions, Equations 9 and 10 give J_E = 0, but J_B < 0. It has been shown by Omura and Nunn (2011) that the component J_B is related to the change of frequency ω′ across one whole subpacket by

(24)

The growth in frequency described by Equation 6 leads to the decrease of S and to the appearance of J_E, which maximizes for S = −S_max. Increase in J_E is followed by growth in amplitude as described by Equation 7. The emission also propagates away from the equator, experiencing further convective growth (Equation 2). The growth in the source is limited by the optimum amplitude (Omura & Nunn, 2011). As was the case with the first chorus equation (Equation 6), we need to include the shift of the source into the definition of the optimum amplitude. Let us introduce the ratio τ = T_N/T_tr of the nonlinear transition time T_N for formation of the nonlinear resonant current, and the nonlinear trapping period

(25)

Now we put forward an assumption that the optimum amplitude for nonlinear growth is reached when the frequency sweep rate over a trapping period ω′/T_N is equal to the sweep rate ∂ω/∂t given by Equation 6. Since S = −S_max in the source, we have J_B = −J_{B, max}J₀, where J_{B, max} ≈ 1.3 can be obtained by numerical evaluation of Equation 10. With this assumption, we can use Equations 25, 24-25, 24, and 6 to obtain the optimum amplitude

(26)

After the wave amplitude reaches B_opt, the nonlinear growth mechanism breaks down. At the same time, the strongest resonant current is released into the upstream. As explained in section 2.2, it forms a helical structure which continually radiates a whistler mode wave at a frequency that matches the frequency of the initial wave at the point where the current has been created, that is, a frequency ω₁ = ω₀ + Δω₁, where ω₀ is the wave frequency of the initial subpacket and Δω₁ is the frequency difference measured at the point where the optimum amplitude was reached (Point 1 in Figure 1). To model a smooth decrease in amplitude of the initial subpacket, we simply switch the sign of the right-hand side of Equation 7. It is further assumed that the new wave, produced by the radiation from the helical current, cannot replace the previous subpacket until its amplitude drops below B_thr (Point 1′′ in Figure 1). Using the group velocity V_g of the whistler mode wave and the resonance velocity V_R of the particles, this corresponds to a wave source located in the distance (Point 1′ in Figure 1)

(27)

starting at time

(28)

The time interval between Points 1′′ and 1′ was denoted Δt₁ = t_end − t_max. Since the radiation emitted by the helical current is coherent, it is immediately subjected to the nonlinear growth effects, provided it reaches the threshold amplitude. A new subpacket is then established at Δh₁ and the process repeats (Points 2, 2′, 2′′ in Figure 1, etc.). The flowchart of our model is sketched in Figure 2.

It will be shown later in section 3.2 that the helical current can indeed be strong enough to emit waves with amplitudes larger than the threshold value B_thr, based on Equations 22 and simulated J_R. The simulation will also confirm that the ratio J_B/B_w from Equation 24 attains large values only near the source, suggesting that the nonlinear frequency growth happens only in that region.

3 Numerical Simulation

3.2 Results

The Equations 1 and 2 are first solved for a set of parameters listed in Table 1 under the row named “Mid”. The chosen value of the magnetic field parameter corresponds to an L shell value of L = 4.5 and equatorial gyrofrequency , where we used the value 3.1 · 10⁻⁵ T for the equatorial strength of the dipole field at the surface of the Earth. The time step is set to and the grid spacing is . In Figure 3 we present time-space plots of the wave frequency ω, wave amplitude B_w, resonant current density J_R and its components J_E, J_B, and the ratio J_B/B_w. According to Equation 24, frequency growth should happen only where the J_B/B_w ratio plotted in Figure 3f is large. This coincides with the source region, supporting thus the validity of our model. Figures 3c and 3d show that while the J_B component of the resonant current density dominates in the downstream, it has values comparable to J_E close to the source region. The peak values of the total resonant current density J_peak in the source range from −0.39 · 10⁻⁴ J_norm (first subpacket) to −1.06 · 10⁻⁴ J_norm (last subpacket), where is a normalization factor. Following the scheme in Figure 1, we take the peak value for the first subpacket and plug it into Equation 22 to calculate the strength of the magnetic field of the newly radiated wave. Assuming the length of one loop of the helix , we get B_tot = 3.2 · 10⁻⁵ B_eq, which we can compare with the local threshold amplitude B_thr = 1.0 · 10⁻⁶ B_eq. The helical current can span over hundreds of loops, seemingly increasing the estimate by up to 2 orders of magnitude. However, due to the frequency growth in the source, the pitch of the helix is changing and so each section radiates at a different frequency, limiting thus the spatial range we can use for our calculations. We will discuss this in more detail in section 5.

Table 1. Table With Input and Output Parameters

		Q	τ
Mid		0.5	0.5	5.0	0.3	0.4	0.15	1.36
Low		0.25	0.25	4.0	0.2	0.3	0.12	0.86
High		1.0	1.0	6.0	0.4	0.5	0.20	3.07
Set 1		0.25	0.25	5.0	0.3	0.4	0.15	1.36
Set 2		0.5	1.0	6.0	0.4	0.4	0.20	0.86	Mid	Set 1	Set 2
NS	Low	13	12	4	7	32	9	31	30	15	67
	High	24	142	30	25	28	29	26
	Low	4,400	1,700	3,300	3,400	6,700	3,700	2,800	3,800	3,700	2,100
	High	1,900	6,500	2,800	2,200	2,500	3,200	6,600
	Low	2.8	13.1	6.8	2.0	5.0	5.3	13.7	7.1	7.4	13.8
	High	12.4	4.8	7.8	11.2	9.8	8.2	2.5
	Low	310	220	30	220	580	100	300	400	400	300
	High	230	590	370	250	300	350	660
	Low	0.6	2.2	0.3	0.4	0.8	0.5	1.5	1.5	1.1	1.4
	High	2.8	0.7	1.3	2.5	1.6	1.6	1.5
	Low	0.290	0.500	0.220	0.247	0.500	0.257	0.500	0.500	0.500	0.500
	High	0.500	0.500	0.500	0.500	0.500	0.500	0.500

• Note. Values in row “Mid” of the upper section of the table were used to produce the results in Figure 3, rows “Low” and “High” show alternate values for each of the parameters, and rows “Set 1” and “Set 2” represent a set of values compiled from the three previous rows. Values in rows “Set 1” and “Set 2” were used to produce the results in Figure 4. The lower section of the table lists values of the following output parameters: number of subpackets N_S, upstream shift of the source h_elm, frequency sweep rate Δω/Δt, the time duration t_elm, the maximum amplitude B_{w, max}, and the maximum frequency ω_max. In this lower section, rows labeled as “Low” (“High”) were obtained from simulations with input parameters from the “Mid” set of input parameters, but in each column we replaced the “Mid” value of the respective input parameter by its “Low” (“High”) value. Values of the output parameters for the three sets of input values “Mid”, “Set 1,” and “Set 2” are shown in the three additional columns on the right side of the table. The sweep rate, the time duration, and the maximum amplitude were always computed at a distance .

To show the effect of the model's parameters on the overall result, we increased or decreased the values of the parameters one by one according to rows “Low” and “High” in Table 1. We recorded the number of subpackets N_S, upstream shift of the source location across the whole chorus element h_elm, the time duration t_elm, frequency sweep rate Δω/Δt, and the maximum amplitude B_{w, max}. Sweep rate, time duration, and maximum amplitudes are calculated for , which is approximately equal to 2,500 km or to a magnetic latitude λ_m = 5° for L = 4.5. If we measured the maximum amplitudes at larger h, they would grow steadily up to unreasonable values (B_{w, max}/B_eq > 0.1), which is caused by the assumption of parallel propagation of whistler modes, which cannot be justified further from the equator, as was shown by systematic analysis of spacecraft measurements (Santolík et al., 2014) as well as by theoretical considerations of chorus propagation in small ducts (Hanzelka & Santolík, 2019).

From a combination of values from the rows “Mid,” “Low,” and “High” in Table 1, two new sets of parameters were assembled, “Set 1” and “Set 2,” with the goal of obtaining a very low and a very high number of subpackets, while keeping the upstream shift, time duration, and maximum wave amplitude of the element at reasonably low values. The first set consists of “Low” values of τ and Q and “Mid” values of the rest of the parameters. The second set consists of a “Low” value of a, “Mid” values of Q and V_⊥0 and “High” values of τ, ω_pe, ω_phe, and U_th, ‖eq. The resulting time-space plots of wave frequencies and amplitudes are presented in Figure 4. With the first set we managed to push the number of subpackets down to N_S = 15, while with the second set a very large value N_S = 66 was obtained.

As we have seen in section 2, most of the simulation parameters influence the model in a highly complex manner. However, with the use of the results presented in Table 1 and Figure 3, we can observe some patterns. The effect of the parameter τ is probably the most obvious, as it is found only in the formula for the optimum amplitude, Equation 26. Low values of τ give large optimum amplitudes, allowing the wave frequency to grow more rapidly within one subpacket, which leads to a lower number of subpackets and that in turn decreases the total upstream shift of the source. The time duration is decreased due to the strong frequency growth as well. And naturally, higher maximum amplitudes in the source result in higher amplitudes in the downstream. The influence of the optimum amplitude on the results is visible also with the altered values of the other model parameters, but it is combined with effects caused mainly by changes in J_E and B_thr. Increase/decrease in Q has the same qualitative effect as equivalent decrease/increase in τ, except for the low number of subpacket for small Q which is caused by the early termination of the simulation due to low values of optimum amplitudes in the upstream. Higher plasma frequency values can significantly decrease h_elm, but they have little effect on the other output parameters. Increased values of the density of hot plasma population, expressed through ω_phe, and perpendicular velocity V_⊥0, affect the results qualitatively in the same way as an increase in Q. Low values of V_⊥0 can strongly increase the drift of the source and the time duration of the element. The parallel thermal velocity has the most complex influence due to its appearance in the exponential in Equation 12 as well as in the denominator of the formula, but the overall trend in the observed resulting parameters is similar to the effect of ω_phe. Finally, magnetic field inhomogeneity parameter a can strongly influence the sweep rate and the drift of the source.

To better understand what the chorus element could look like in the measurements of a stationary spacecraft, we plot the time evolution of wave frequency and amplitude in Figure 5 for the three sets of parameters “Mid,” “Set 1,” and “Set 2”. The position in space is fixed to latitudes of 5° (red lines) and 15° (blue lines). In Figures 5a and 5c we can clearly see frequency drops between adjacent subpackets, while in panel (e) this behavior becomes indistinct due to the large number of subpackets in the fine structure. Also, with rising frequency the subpackets are getting shorter and the ratios between the increase and the following drop in frequency within one subpacket are decreasing. The envelope of the amplitudes follows the dependence of the optimum amplitude on frequency (see, e.g., Omura & Nunn, 2011, Figure 3a for comparison). With rising frequencies the peaks in the amplitude plot are getting smoother due to increasing dispersion of the whistler mode waves propagating in cold plasma. Dispersion also causes decrease of the relative height of the peak (from base to top), making the fine structure more homogeneous.

Last but not least, we have tested the influence of the initial value of the wave amplitude of each subpacket. We determined that as long as the threshold amplitude B_thr is by at least 1 order of magnitude smaller than the optimum amplitude B_opt, any initial amplitude that ranges from about 1.5 B_thr to 3.0 B_thr has negligible effect on the results of the simulation. Similarly, decreasing integration steps in space and time by half did not lead to any changes in the values of output parameters.

4 Comparison with Observation

High quality electromagnetic wave measurements provided by the two Van Allen Probes were used to identify large-amplitude chorus events in the radiation belts. One such event, detected by the Van Allen Probe B spacecraft on 12 September 2014, is presented in Figure 6. Figures 6a and 6b respectively show the frequency-time power spectrograms obtained from the magnetic field and the electric field measurements, recorded by the EMFISIS Waves instrument (Kletzing et al., 2013) in the morning sector at McIlwain's L = 5.61 and magnetic latitude λ_m = 5.24° northward from the magnetic equator. A sequence of intense chorus elements is clearly seen in both spectrograms below one half of the local electron cyclotron frequency, which is shown as a white or black solid line on the spectrograms. These electromagnetic waves have a right-hand circular polarization, indicating the presence of the whistler mode in Figure 6c obtained using the method of Santolík et al. (2002).

The planarity of the magnetic polarization obtained by the singular value decomposition (SVD) method (Santolík et al., 2003), plotted in Figure 6d, is above 0.8 in the lower-frequency parts of the elements between 1.2 and 1.5 kHz, corresponding to the presence of a single plane wave in a given frequency-time bin of the spectrogram. The planarity is below 0.8 in the upper frequency parts of the elements extending up to a frequency of 1.7 kHz, suggesting that the plane wave approximation should not be used above 1.5 kHz.

The angle θ_k between the wave vector and local magnetic field line is lower than 10–20° below 1.5 kHz, as shown in Figure 6e. The higher values observed at larger frequencies are not reliable under the plane wave assumption. The azimuth of the wave vector in Figure 6f shows a predominant outward propagation in the plane of the local magnetic meridian. Finally, Figure 6g shows that spectral estimates of the Poynting vector, obtained using a method of Santolík et al. (2010), give directions outward from the magnetic equator.

The data recorded in the burst mode of the EMFISIS Waves instrument have a sampling rate of 35 kHz and a 16-bit dynamic range, allowing thus for a detailed analysis of the fine structure of chorus. Figures 7a and 7b show detailed waveforms of the first chorus element from Figure 6. The analysis method used in Figures 7c–7e is similar to the method used for measurements of the Cluster mission by Santolík et al. (2004) and the same as the analysis procedure used for another interval of Van Allen Probes measurements by Santolík et al. (2014): The calibrated waveform is pass-band filtered between 0.4 and 3 kHz, and analytic signals are constructed using the Hilbert transform. Their instantaneous amplitudes are shown in Figure 7c. The instantaneous frequencies plotted in Figure 7d are obtained as time derivatives of the phases of the complex analytic signals, while both the instantaneous phases and amplitudes are used to construct instantaneous spectral matrices, whose SVD analysis provides us with estimates of the instantaneous wave vector angles plotted in Figure 7e.

The analyzed chorus element is composed of subpackets, in consistence with the assumptions made in the model described in section 2.3. The instantaneous frequency is globally rising with time but sometimes it steps back at the boundaries of the subpackets. This is consistent with the simulation results in section 3. The input and output parameters analyzed in Table 1 cannot be readily compared with the observation since we do not measure Q and τ, which have both strong influence on the output parameters. Also, the assumption of bi-Maxwellian distribution, included in Equations 12 and 13, need not hold, making the parameters V_⊥0 and U_th, ‖eq hard to interpret. Nevertheless, we can still look at the properties of the analyzed element and see that the parameters N_S ≈ 23, Δω/Δt ≈ 1.8 kHz/s, and t_elm ≈ 400 ms are within a multiple of 2 from the output parameters obtained in the simulation with a = 3.07 · 10⁻⁷ c⁻²Ω_e0, which corresponds to L = 5.5.

Figure 7 clearly shows that the waveforms of the perpendicular and parallel components behave differently, their subpacket structure is different and their estimated instantaneous frequencies are also slightly different. This is strongly reflected by the instantaneous wave vector angle which changes its value within each subpacket. As it was already noted for another case analyzed by Santolík et al. (2014), the amplitude maxima generally correspond to the minima of the instantaneous wave vector angle.

5 Discussion and Conclusion

In the development of our model of the fine structure of rising tone chorus emission, we decided to base it on the nonlinear growth theory described in Omura et al. (2008) and the follow-up papers. There exists another prominent theory of the chorus emission, summarized, for example, in Trakhtengerts (1999), which is based on the backward oscillator regime of cyclotron masers in space. It has been successfully applied to explain the time delay between chorus elements and their frequency sweep rate (Demekhov, 2017, and references therein), but it has not yet explained their fine structure.

A crucial role in the subpacket formation process is played by the electromagnetic radiation emitted from the resonant electrons leaving the source region. We have shown that the emitted wave should be theoretically far above the threshold amplitude, possibly even reaching the optimum amplitude, which would stop the nonlinear growth mechanism. However, the computation relied on the current having a shape of a perfect helix. In reality, the magnitude of the current is dependent on the gyrophase bunching process. Without gyrophase bunching of the resonant electrons, there is no net current. Also, the whole process of the release of the current at h_i and radiation of a new wave at h_i + 1 is a discretization of a continuous motion of a radiation source with a spatial spread, which we could not describe due to simplifications that appear in the nonlinear growth theory, namely the point-like source approximation. Future studies should focus on the calculation of the magnitude of the resonant current. Natural extension of results of our simplified wave field model would be to use it for analysis of the velocity distribution of test particles propagating through the wave field modeled in this paper. In full-particle simulations, the radiation appears naturally in the solution of Maxwell equations for the particle system.

Another effect that can decrease the power of the emitted wave is the changing pitch of the helix. As the frequency of the wave inside the primary subpacket continuously increases, the helical current must copy the structure and change its pitch. This would lead to broadening of the spectral peak of the emitted wave, and to decrease of its maximum power. Since the amplitude of the current in the source has a peak (see Figure 3e, also compare with amplitude peaks in Figure 5b which partially copy the evolution of current), we do not expect this effect to be very prominent. Nevertheless, it is clear that the true nature of this radiation process is more complex than shown in our model. A more consistent approach to the antenna effect can be found in Trakhtengerts et al. (2003), where they compute the radiated power and frequency shift directly from the transport equations for the wave amplitude and nonlinear phase. Since they do not consider any subpacket structure, the antenna length becomes much longer and dephasing starts to play a major role. They conclude that the frequency shift due to the antenna effect should be about 100 Hz in typical magnetospheric conditions, which is similar to our result, and that the amplitude of the new wave B_w/B_eq is between 10⁻⁵ and 10⁻⁴, which is above the threshold amplitudes considered in this paper.

The comparison of simulation results with observations of the Van Allen Probes confirms that the drops in frequency between subpackets, which are a fundamental part of our model, can be observed in large amplitude chorus elements. The upstream shift of the source, which is another important feature of the model, cannot be determined from measurements of a single spacecraft, but indirect indications of a similar effect have been reported by Taubenschuss et al. (2017) for bidirectional chorus wave packets. Two satellites with a small spatial separation (hundreds of kilometers) should be in principle able to directly intercept one chorus element inside the source at different stages of its development. If this proposed drift of the source were real, one satellite (at the equator) would see the whole frequency range of the element, while the other one (shifted slightly upstream) would see only the upper part of the range, and the first coherent, large-amplitude emission would appear with a significant time delay with respect to the first satellite's measurement. Short distances between spacecraft with highly sensitive wave instruments were achieved during several close separation campaigns of the four-spacecraft Cluster mission (see, e.g., Němec et al., 2014), and additional work is needed to identify signatures of this effect for special configurations when different spacecraft are located close to a single magnetic field line, at transverse separations lower than a typical transverse size of generation regions of separate chorus wave packets, that is, on the order of 100 km according to Santolík and Gurnett (2003) and Santolík et al. (2004).

The only simulation that clearly showed and analyzed a shift of the source region within a nonlinear theory was the simulation of EMIC waves by Shoji and Omura (2013), where the upstream drift of the source was qualitatively similar to our chorus simulation, but we cannot make any quantitative comparison due to the different nature of the whistler waves and ion cyclotron waves. Some less well-behaved movement of the source has been observed in chorus simulations as well, for example, in the full-particle simulations of Hikishima and Omura (2012), but it was not properly discussed there.

Another point that must be mentioned in the discussion of our results is the choice of ranges of parameter values which we used in simulations. While the field inhomogeneity a is given by the dipole field model and plasma frequency ω_pe can be chosen based on measurements in the equatorial region of the outer radiation belt, the choice of the remaining parameters is less obvious. The most important constraint imposed on the parameters is that B_thr ≪ B_opt must hold for the initial frequency. Our goal was to keep the values of ω_phe, V_⊥0, and U_th, ‖eq as low as possible, because in general, very hot and dense distributions are less likely to occur. Since all our simulations started at frequency ω = 0.2 Ω_e0, that is, at a fairly low value, we had to settle for hot plasma frequency of about 0.3 Ω_e0, which corresponds to relative density n_hot/n_cold = 3.6 · 10⁻³ for ω_pe = 5.0 Ω_e0. This is because the ratio B_thr/B_opt increases rapidly as the wave frequency decreases, as was shown by Omura and Nunn (2011). Even with these high hot electron densities, a small change of parameters could lead to large drifts of the source, which can cause the optimum amplitude to decrease below the threshold amplitude. This is demonstrated in Table 1, where the maximum frequency ω_max does not always reach the limiting frequency ω_fin, resulting in very short chorus elements. The quantities Q and τ are essentially free parameters of the nonlinear growth theory, since they cannot be estimated without performing a self-consistent simulation, and therefore can be used to tweak the results to certain extent.

One of the consequences of the rather high values of hot plasma density are the large overall amplitudes of resulting whistler waves, reaching typically a few percent of the background magnetic field (Figure 5). These results are overestimated because we have limited our study to parallel propagation. The θ_k values can also reach tens of degrees inside the source region (Santolík et al., 2009). Even in cases where the propagation is globally quasi-parallel (Figure 6) the θ_k values vary at time scales of subpackets (Santolík et al., 2014), as we can also see in Figure 7. Energy of oblique whistler waves is transferred back to electrons through the Landau resonance (Hsieh & Omura, 2018), decreasing thus the observed wave amplitudes. The two-dimensional nature of the chorus emission also has significant influence on the particle acceleration, as was shown by Omura et al. (2019). Simulation of the two-dimensional evolution of whistler emissions is a hot topic in space plasma science—PIC simulation in curvilinear coordinates by Ke et al. (2017) show changes of wave normal angle inside the subpackets, but the whistler waves are still quasiparallel. The 2-D PIC simulations of Ke et al. (2018) show that the parametric decay of a parallel whistler wave can result in creation of daughter whistler waves with high obliquity, which could possibly explain the wave normal angle values observed in Figure 7e. Crabtree et al. (2017) even suggest that the chorus generation mechanism is inherently three dimensional, as they discovered a smooth change in the azimuthal angle of the wave vector within single subpackets.

To summarize, we have shown that a model based on the nonlinear growth theory and the antenna effect can be used to simulate growth and propagation of single chorus elements with subpacket structure. The model features steep drops in frequency at the point where one subpacket transitions to the next one and an upstream drift of the source region with increasing wave frequency. The first feature was confirmed by observations of the Van Allen Probes spacecraft, the second one appears in self-consistent particle simulations. Time duration and frequency sweep rate of the element and the number of subpackets obtained through simulations are comparable to those observed in a typical event of intense chorus recorded by the Van Allen Probe B spacecraft. The model can be used in test particle simulations (see Nunn & Omura, 2015) to determine the effect of subpackets on particle acceleration—this is left for future studies.