Whistler Mode Waves Below Lower Hybrid Resonance Frequency: Generation and Spectral Features

Equatorial noise in the frequency range below the lower hybrid resonance frequency, whose structure is shaped by high proton cyclotron harmonics, has been observed by the Cluster spacecraft. We develop a model of this wave phenomenon which assumes (as, in general, has been suggested long ago) that the observed spectrum is excited due to loss cone instability of energetic ions in the equatorial region of the magnetosphere. The wavefield is represented as a sum of constant frequency wave packets which cross a number of cyclotron resonances while propagating in a highly oblique mode along quite specific trajectories. The growth (damping) rate of these wave packets varies both in sign and magnitude along the raypath, making the wave net amplification, but not the growth rate, the main characteristic of the wave generation process. The growth rates and the wave amplitudes along the ray paths, determined by the equations of geometrical optics, have been calculated for a 3‐D set of wave packets with various frequencies, initial L shells, and initial wave normal angles at the equator. It is shown that the dynamical spectrum resulting from the proposed model qualitatively matches observations.


Introduction
This paper presents a theoretical analysis of wave phenomena, an example of which is shown in Figure 1. The figure displays the overview time-frequency power spectrogram of magnetic field fluctuations recorded by the STAFF instrument on board the Cluster 3 satellite. The frequency range of the emission falls in the interval between proton gyrofrequency and lower hybrid resonance (LHR) frequency; it has been observed in the near-equatorial region on L shells from about 4 to 4.4. Emissions in this frequency-space domain, referred to as equatorial noise, have been known for years (Russell et al., 1970). Basic ideas relating this wave phenomenon to the interaction at high proton cyclotron harmonics have been put forward by Gurnett (1976). Observations of equatorial noise in conjunction with warm trapped plasma based on the measurements from Dynamics Explorer 1 (DE 1) have been presented by Olsen et al. (1987). More details on this matter may be found in Santolík et al. (2002).
As is well known (Ginzburg & Rukhadze, 1972), in a cold plasma, there is only one propagating mode in the indicated frequency band, which is generally called "whistler mode." In the electron-proton plasma, this mode is a continuation of magnetosonic wave mode from frequencies below the proton cyclotron frequency to higher frequencies, through the LHR frequency up to the min( c , p ), where c , p are the electron cyclotron and plasma frequencies, respectively. The equatorial noise is observed in the magnetospheric region where protons are dominant, but not the only ions. Although in this case the whistler mode is decoupled from the magnetosonic mode, these whistler mode waves that form the equatorial noise are often also referred to as magnetosonic waves (abbreviated as MSW or MS waves further on). In what follows, we will use either of these terms. Concerning the terminology, two remarks are in order. First, in multi-ion species plasma the magnetosonic branch, which starts below cyclotron frequency of the heaviest ion, has a resonance (refractive index tends to infinity) near the cyclotron frequency of the second heaviest ion, so, strictly speaking, the low-frequency whistler mode waves below the LHR frequency cannot be called magnetosonic waves anymore. And second, the term equatorial noise is used to describe the wave phenomenon, but not the plasma mode, and is independent of how the wave mode is called.
It has been suggested by several authors that equatorial noise is generated by an unstable ring-like distribution of energetic protons (see, for example, Perraut et al., 1982). Simultaneous observations of proton ring Figure 1. The dynamical spectrum of the B z magnetic field component as observed by STAFF instrument on board the Cluster 3 spacecraft on the 6 of July 2013. The 15th, 20th, 25th, and 30th harmonics of the proton cyclotron frequency are indicated by solid black lines. The figure displays the same event as has been investigated by Balikhin et al. (2015) and Walker et al. (2015) (see Figure 2 from Balikhin et al. (2015) and Figure 2 from Walker et al. (2015)). distribution and MS waves by Van Allen Probes during small storm on 14-15 April 2014 have been reported by Xiao et al. (2015). However, previous considerations suffered two disadvantages: they used oversimplified distributions and were limited to the analysis of the wave growth rates. In the magnetosphere, where the wave normal vector and plasma parameters vary drastically along the wave packet trajectory, the value of the growth rate calculated for fixed wave normal vector and fixed ratio of the wave frequency to the proton gyrofrequency makes little sense, and it is necessary to calculate the wave packet energy density, proportional to the wave field amplitude squared, with the account of parameter variations along the ray trajectory.
In this paper we follow up the idea about generation of the equatorial noise by an unstable proton distribution but perform calculations of the wave energy density for fixed frequency wave packets and, thus, electromagnetic field spectral density, taking into account the variation of the wave normal vector and plasma parameters along the ray trajectory, using a realistic distribution function of energetic protons. In calculating the proton contribution to the growth rate of magnetosonic waves, we use the approach that had been put forward more than 50 years ago (Dawson, 1961;O'Neil, 1965) and has since then been used in a great number of works on wave-particle interactions in the magnetosphere (see, e.g., Karpman et al., 1975;Meerson et al., 1979;Trakhtengerts & Rycroft, 2008, and references therein). This approach consists of the following. The magnetospheric plasma is assumed to be composed of two parts: a "cold" bulk of electrons and ions that determines the plasma dispersion relation and low-density suprathermal electrons and ions which participate in resonant interaction with the wave and are responsible for wave growth or damping. Apart from this, if the wave growth (or damping) rate is less than the inverse time of resonant interaction, determined by plasma inhomogeneity, which we will also assume, the resonant particle distribution function can be found in adiabatic approximation with respect to the wave amplitude, that is, neglecting the amplitude variation during the time of resonant interaction (O'Neil, 1965) (see also Shklyar and Matsumoto (2009), and references therein).
The next section offers some general considerations of MSW propagation and generation. In section 3, the dispersion relation and the growth rate of these waves are discussed. Section 4 deals with the calculation of energy density along the wave packet trajectory and with spectrum modeling. The results of the work are briefly summarized in section 5. We should mention that the waves under discussion, which are generated by an unstable proton distribution, are believed to play an important role in the dynamics of energetic electrons of the radiation belts (see, e.g., Horne et al., 2007;Li et al., 2014;Mourenas et al., 2013, and references therein). This problem, however, is beyond the scope of the present study which is devoted exclusively to generation and spectral features of these waves.

General Remarks on MSW Propagation and Generation
To understand the peculiarities of the magnetosonic wave spectrum observed on Cluster, the following points are of importance.
When discussing the waves with frequencies less than the LHR frequency, in particular, in the range of high proton cyclotron harmonics, both electron and ion contributions to the wave dispersion, as well as the wave growth and damping should be taken into account. For electrons, the Cerenkov and the first cyclotron resonances are the most important, while all other cyclotron resonances correspond to much higher resonance velocities, thus much smaller amount of resonant particles. Among those two, the Cerenkov resonance, for which the resonance velocity is v r0 = ∕k ∥ , usually leads to the wave damping. Here and further, k ∥ denotes the parallel component of the wave normal vector. The first cyclotron resonance, corresponding to the resonance velocity v r1 = ( − c )∕k ∥ , may lead to wave growth, provided that plasma is unstable due to the presence of a loss cone or/and temperature anisotropy in the electron distribution. However, for the small frequencies under discussion, that is, ≪ c , |v r0 | ≪ |v r1 |; thus, these resonances cannot compete, and the waves are not excited by the electrons. Moreover, even if the wave is excited by protons, it will be damped by electrons, unless k ∥ is so small that |v r0 | ≫ v Te , where v Te is a thermal velocity of energetic electrons. Thus, the first conclusion is that only waves, having the wave normal angle close to ∕2, may be observed. Such waves can be excited by an unstable proton distribution at high-order cyclotron resonances in the regions where the wave frequency and the proton cyclotron frequency Ω satisfy the relation ≃ nΩ for some integer n, since in this case the resonance velocity corresponding to the nth proton cyclotron resonance, that is, V Rn = ( − nΩ)∕k ∥ remains finite. We should stress that the resonant interaction between low-frequency whistler mode wave and protons at high-order cyclotron resonances is efficient only under the condition that proton transversal velocity V ⟂ is large enough, namely, V ⟂ > ∕k ⟂ , where k ⟂ is a perpendicular component of the wave normal angle (Shklyar, 1986); thus, an efficient resonant interaction between the wave and thermal protons, which could potentially result in heavy wave damping, does not take place even if the resonance velocity tends to zero.
Another important point concerning the generation of magnetosonic waves via high-order proton cyclotron resonances consists of the following. As all low-frequency whistler mode waves, these waves propagate (in the sense of their group, but not phase velocity) almost along the ambient magnetic field (Storey, 1953), except close vicinities of the reflection points where k ∥ = 0 and, thus, V Rn → ∞. Then, in order that a wave remained in the nth cyclotron resonance with energetic protons for a sufficient time and experienced a significant amplification, the proton cyclotron frequency should change slowly along the wave trajectory, which is possible only in the equatorial region. Thus, only the waves propagating in the equatorial region at high wave normal angles can be observed in the frequency band under discussion. A typical trajectory of such magnetosonic wave is shown in Figure 2. The parameters of wave propagation along this trajectory that are essential for our discussion are illustrated by Figure 3. For waves starting with the wave normal vector directed from the Earth, the L shell increases along the ray trajectory, while the latitude oscillates around the equatorial value = 0. The index of refraction increases monotonically, while the wave normal angle oscillates around ∕2. Parallel group velocity V g∥ oscillates around zero in accordance with the wave packet bounce oscillations around the equator, while perpendicular group velocity V g⟂ decreases in magnitude remaining of the same sign.
In the present study, we limit ourselves to wave propagation in meridional plane. This permits to find parallel and perpendicular components of the wave normal vector that are essential for calculating the wave amplification, which is the main focus of our work. A 3-D propagation of MS waves in relation to their MLT distribution has been studied by Xiao et al. (2012).
Variations of the quantity ∕Ω, the integer part of which determines the order of the closest proton cyclotron resonance, as a function of time along the trajectories of quasi-monochromatic wave packets, are shown in Figure 4. We see that at least on the average, the cyclotron frequency decreases along the ray trajectory and, thus, the order of cyclotron resonance crossed by the wave packet increases. One important consequence of this observation, related to the lower frequency part of the spectrum where the decrease of cyclotron frequency along the ray trajectory is monotonic, consists of the following. Wherever the wave is observed, it came from the region where the proton gyrofrequency was larger than the local gyrofrequency. It means that if at the observation point the wave frequency is slightly above the nth cyclotron harmonic, the wave should have come from the region where its frequency was equal to nΩ. As we will see below, in this region the quantity that characterizes the wave growth ( > 0) or damping ( < 0) is always negative. On the contrary, if at the observation point the wave frequency is below the nth cyclotron harmonic, the corresponding wave packet did not cross the exact proton cyclotron resonance right before the observation. Thus, one may expect the spectral intensity below the local proton cyclotron harmonic to be higher than above the corresponding frequency. This spectral feature seems to be consistent with the observations. We should stress that this statement is by no means categorical, as the wave packet amplification is gathered along the whole raypath, but not only in close vicinity of the observation point; and the exact cyclotron resonance, although making the main contribution, is not the only one that contributes to the wave growth (damping) rate.

Model of the MSW Generation
We will assume that MSW spectrum is excited due to the loss cone instability. The approach that we will use to describe the excited emission consists of the following. We will suppose that the wave field can be presented as a sum of constant frequency wave packets propagating according to the equations of geometrical optics, with the corresponding variation of the wave normal vector along the ray trajectory, and the amplitude growth or damping determined by the local linear growth (damping) rate L , calculated with the local values of the wave normal vector and plasma parameters. In order to implement this approach, we need to define the wave dispersion and polarization, and the wave growth (or damping) rate based on a certain unperturbed distribution function. We will give the corresponding relations without detailed derivation, as they have been presented in a number of papers (see, e.g., Shklyar & Matsumoto, 2009;Shklyar et al., 2004;Walker et al., 2015).

Dispersion Relation and Wave Polarization
The dispersion relation for whistler mode waves which remains valid below the LHR frequency has the form where LH is the LHR frequency, and k ∥ = kcos , k ⟂ = ksin , where is the wave normal angle, that is, the angle between the wave normal vector and the ambient magnetic field. The equations of geometrical optics used in this study are based on the dispersion relation (1) and have the usual form where (k, r) is determined by (1). The details about their solution may be found in Shklyar and Jiȓícek (2000) and Shklyar et al. (2004).
In this study, we use the dipolar model of the ambient magnetic field, and gyrotropic model of cold plasma density in which 2 p ∝ c , with p ∕ c = 6 on L = 3 at the equator. As the calculations show, including a plasmapause, or/and another model of plasma density, for example, 2 p ∝ c , ≠ 1 does not change the results qualitatively.
The electric and magnetic fields of the wave packet are assumed to have the form where the local wave normal vector is given by k = ∇ , E is a complex scalar amplitude of the wave electric field, and a, b are complex polarization vectors of the wave electric and magnetic field, respectively, related by the Farday's law In the local coordinate system with the ambient magnetic field directed along the z axis, and the wave normal vector lying in the (x, z) plane, the components of the polarization vectors a, b have the form where 1 , 2 , and 3 are components of the dielectric tensor k x , k z are the corresponding components of the wave normal vector, and k 2 = k 2 x + k 2 z . In the frequency band of interest, components of the dielectric tensor 1,2,3 can be written in the form where pi is proton plasma frequency and UH = √ 2 p + 2 c is the upper hybrid frequency. From (7) it follows that for frequencies under consideration 1 < 0 and 2 < 0, thus, for the wave normal angles close to ∕2, a z and a x satisfy the relation a z ≪ a x . Apart from this, the following relations hold: | 1 | ≪ | 2 | ≪ | 3 |, so that the wave electric field is almost linearly polarized along the x axis, although it has a small y component corresponding to right-hand polarization. As for the polarization of the wave magnetic field, from (5) it follows that b x and b y components are close to zero, so b z component is the largest one. It means that the magnetic field of the waves under discussion is almost linearly polarized along the ambient magnetic field.

Unperturbed Distribution Function of Energetic Protons
Turning to the linear growth rate due to interaction with energetic protons, we first specify their unperturbed distribution function. This function should depend on the constants of the motion in the absence of the wave field, for which we will choose particle energy W and equatorial pitch angle determined by where = MV 2 ⟂ ∕2Ω is proton magnetic momentum and Ω 0 is the equatorial proton gyrofrequency. We will use a simple distribution of the form: where C is the normalization constant defined below. The distribution function (9) depends on two parameters: the loss cone parameter > 0 and the effective temperature T. Being considered as the function of W and , the distribution (9), for any fixed W, tends to zero when → 0 and has maximum at = ∕2. The normalization constant C can be expressed through the density of energetic protons, n h . Integrating (9) over velocity space we obtain where The function I( ) obtained by numerical calculation of the integral (11) is shown in Figure 5. For the chosen form of distribution function, the density of the energetic protons is not constant, varying as (Ω 0 ∕Ω) (see (10)). Along with n h , we will further use the constant quantity n h0 = C v 3 T I( ) equal to the energetic proton density at the equator.

Linear Growth Rate
For the waves under discussion, the linear growth (damping) rate due to resonant interaction with energetic protons is given by (Shklyar & Matsumoto, 2009 where summation is over all cyclotron resonances, U ∝ |E| 2 is the wave energy density, f 0 is an unperturbed proton distribution function expressed through the particle kinetic energy W and magnetic momentum , and

10.1002/2017JA024416
where V rn is the nth cyclotron resonance velocity for protons: The wave energy density U expressed through the amplitude of the x component of the wave electric field may be written as The quantity u can be found from the general expression for wave energy density (Shafranov, 1967) using the expressions (7) and (4). The result reads We will also give an explicit relation between a x and b z , which allows to express the wave energy density in terms of B z component of the wave magnetic field The quantity A n which enters the expression for the growth rate is the effective amplitude of resonant interaction at the nth cyclotron resonance given by where J n ( ) and J ′ n ( ) are, respectively, the Bessel function and its derivative with respect to the argument . Calculating the derivative (13) with the distribution function (9) and substituting the result into (12), we obtain where Given the initial conditions for the wave packet, that is, the initial coordinate r 0 and initial wave normal vector k 0 , the equations of geometrical optics together with equations (18) and (19) give the current values r, k, and the growth rate along the trajectory of the wave packet.
Proceeding to an analysis of the expression (19) we first notice that the contribution to the growth rate from the nth cyclotron resonance is essentially determined by the exponential factor before the integral. Using wave and plasma parameters typical of the case under discussion (see Figure 3), namely, k ∥ ∼ 1.5 ⋅ 10 −7 cm −1 , and Ω ∼ 32 rad/s and T ∼ 35 keV, we find which shows that only the resonances with n close to the integer part of ∕Ω can contribute to the growth rate. In numerical calculations we used the three closest resonances. It is easy to see that at the point where = nΩ, the nth cyclotron resonance always makes a negative contribution to . Indeed, at the exact nth cyclotron resonance, the expression in square brackets in (19) is equal to − ∕T and, thus, is negative for all , while all other factors are positively defined. On the contrary, for finite V rn and small , the expression in square brackets is positive, and so may be the integral, thus, the growth rate .  4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29  For fixed plasma parameters, that is, at a fixed coordinate r, and the chosen type of unperturbed distribution function (9) and the dispersion relation (1), the growth rate is the function of the wave frequency and the wave normal angle; thus, to show this function as a usual 2-D plot, one of these parameters should be fixed. The dependence of the normalized growth rate ∕Ω on the normalized frequency ∕Ω for the wave normal angle = 89 ∘ is shown in Figure 6. In agreement with the analytical consideration, at exact cyclotron resonances is always negative. Figure 6 corresponds to the loss cone parameter = 1.5. When increases from 1 to 2.5, while all other parameters remain unchanged, the general shape of the increment as the function of wave frequency remains the same, while the magnitude of slightly increases. The dependence of the growth rate on the wave normal angle for four frequencies is shown in Figure 7. The growth rate as the function of frequency and the wave normal angle is presented in Figure 8 by color code. We see that the growth rate depends crucially on both the relative frequency ∕Ω and the wave normal angle , while both of these quantities vary along the raypath. Thus, the local growth rate itself does not characterize the wave packet amplification, which should be considered on its own account.
The growth rate obtained in the present study differs from the growth rates which are encountered in literature (see, e.g., Ma et al., 2014;Zhou et al., 2014) and have the shape of a comb with sharp maxima at exact proton cyclotron harmonics and equal to zero in between. Our growth rate is negative at exact cyclotron harmonics, that is, when = nΩ (cf. Figure 6), a negative contribution to from the nth term in (19) being proved analytically. The reason for the discrepancy might be the following. The expression for the growth rate from the paper by Chen et al. (2010), to which the papers cited above refer, contains |k ∥ | in denominator, similar to expressions (12) and (19) that we use. In the case of interest, that is, for k ∥ → 0, one should be very careful with such an expression. As was mentioned Local growth rates as a function of wave normal angle for four frequencies. Also indicated are the quantities ∕Ω that define the order of the closest cyclotron resonance. All parameters are the same as in Figure 6. above, the magnitude of the growth rate is determined to a considerable degree by exponent For k ∥ → 0 and ( − nΩ) ≠ 0 this factor is exponentially small, and tends to zero despite k ∥ in the denominator of preexponential factor. On the other hand, if k ∥ is small, but finite, and ( − nΩ) → 0, then is a large (in magnitude) negative quantity, at least for the distribution function that we consider. And if both quantities tend to zero simultaneously, there appears indeterminacy of the type 0∕0 in the exponent, and the additional zero in denominator, so that the result depends on the way these quantities tend to zero. However, the question of the correct evaluation of the growth rate in this case is not a matter of a correct limiting process. Neither the expression given by Shklyar and Matsumoto (2009) nor by Chen et al. (2010) is applicable in the case k ∥ → 0, while in the case of ( − nΩ) = 0 and finite k ∥ both expressions remain valid. A different sign of at = nΩ and k ∥ → 0 is probably due to different distribution functions, although another distribution function that we checked, with empty loss cone and logarithmic growth beyond the loss cone boundary also gives negative in this case. Anyway, we believe that basing the consideration on the growth rate calculated for k ∥ → 0 is not pertinent, as this quantity varies along the wave packet path. We should mention that the particular case k ∥ = 0 for homogeneous plasma in which k does not change has been considered by Shklyar (2009).
The linear kinetic dispersion theory, including calculation of the growth rate, for lower proton cyclotron harmonics has been developed by Gary et al. (2010) using subtracted Maxwellian model for nonthermal proton distribution. The authors have shown that the growth rate has maxima at frequencies close, but not exactly equal, to proton cyclotron harmonics and calculated the growth rate as a function of wave normal vector k, for various wave normal angles that provide maximal growth rate close to the corresponding cyclotron harmonic. Since Gary et al. (2010) considered the case of homogeneous plasma and did not calculate the amplification coefficient, they displayed only positive values of the growth rate , but there is no doubt that for the values of k (and, thus, , as the wave normal angle has been fixed) for which is not shown, the latter is negative. Although the growth rate of Gary et al. (2010) has been calculated as a function of k for various wave normal angles, while we calculate the growth rate as a function of frequency for fixed wave normal angle, the qualitative agreement between our results is beyond any doubt. In particular, in both cases, the maximal growth rate decreases with increasing harmonic number, while the frequency range of positive increases. Thus, our results may be thought as an extension of the results by Gary et al. (2010) to higher cyclotron harmonics, up to the frequencies of the order of LHR frequency that are of interest in relation to equatorial noise. Apart from extended frequency range, we also use another type of unstable proton distribution function.

Wave Packet Energy and MSW Spectrum
Modeling of MSW spectrum presented in this section is based on the following assumptions. Ion distribution is unstable in the region L > 3, where the exact figure is taken for the sake of definiteness. The problem is stationary, while the evolution of spectrum along the satellite orbit is due to spatial spectrum variation. The observed spectrum is formed by constant frequency wave packets that propagate from lower L shells to the observation point and are amplified in the unstable region. The observed spectral intensity at a given frequency f is proportional to the energy density of the wave packet with the same frequency.
In a stationary case, the equation for the wave packet energy density U has the form where is the growth rate and v g is the group velocity. Equation (20) may be written in the form where s is the length along the ray considered, is the cross section of a thin ray tube centered on this ray, and v g is the magnitude of the group velocity. Introducing the quantity that is, the energy flux across the ray tube, and taking into account that according to equations of geometrical optics, ds∕v g = d , where is the time of wave packet propagation along the given ray, we rewrite equation (21) as An important realization concerning equation (23) consists in that its formal integration in the case of sign varying would lead to a physically incorrect result because, in linear approximation in use, the wave growth starts from some finite thermal level; thus, the wave energy density cannot decrease below this level, while equation (23) permits such a decrease. That is why in numerical calculations instead of (23) we solved the equation Obviously, all wave packets that form the observed spectrum start in an unstable region of the magnetosphere, that is, at L > 3 under our assumption. The simplification that we adopt consists in that all rays start at the equator from the interval [3, 3 + ΔL(f )] of initial L shells. The quantity 3 + ΔL(f ) is equal to the L value at which the ray starting at L = 3 at the equator crosses the equator again for the first time. As the ray tracing calculations show, this quantity depends on the wave frequency, but it depends very little on the initial wave normal angle. The quantity ΔL is equal to [0.82, 0.47, 0.335, 0.26, 0.205] for frequencies [50,100,150,200,250] Hz, respectively, and we use cubic interpolation for other frequencies. We have performed numerical simulations for 81 linearly spaced frequencies in the interval 50-250 Hz. The waves start from the interval ΔL(f ) with the step L = 0.1. For each frequency and initial L shell we consider six rays with initial wave normal angles 0 = [87.5 ∘ , 88.5 ∘ , 89.5 ∘ , 90.5 ∘ , 91.5 ∘ , 92.5 ∘ ]. For all combinations of frequency, and each initial parameters indicated above, we solve the equation (24) along the ray trajectory, from the starting point up to the satellite orbit corresponding to observations shown in Figure 1, namely, 0.3162L + tan − 1.3093 = 0, where is geomagnetic latitude in radians. The value of P at the satellite orbit corresponding to a given frequency f is associated with the observed spectral intensity at the same frequency.
As one can see from Figure 1, both L and vary monotonically along the satellite orbit, approximately according to the above given relation. We parametrize the satellite orbit by L value and choose the width of the grid over L equal to ( L) grid = 0.01, while the coordinate L varies from 3.92 to 4.39 on the satellite orbit under discussion. We choose L as one coordinate of the net on which the spectral intensity is calculated, the other one being the frequency f , of course. The (L, f ) grid contains 3,807 cells. If the ray trajectory corresponding to the wave packet of the frequency f reaches the satellite orbit at some point (L, ), then the calculated value of P falls into corresponding (L, ) sell. The number of ray trajectories that we calculate is not given by simple multiplication since equidistantly divided initial L span depends on frequency. In fact, the model described above deals with 2,142 rays. Since the number of rays is less than the number of cells, and since not all rays reach the satellite orbit, in most cases each cell contains only one value of P or is empty. However, if more than one value of P falls into the same cell, they are summed up. The spectrogram of P values obtained in the way described above is shown in Figure 9. We should stress that the quantity P(f ) is proportional, but not equal, to the power spectral density of the wave magnetic field displayed in Figure 1. Simulation of the power spectrogram of the wave magnetic field, which requires the evaluation of the ray tube cross section along wave packet trajectories is beyond the scope of the present work.

Summary
Cluster measurements of magnetic field fluctuations below the LHR frequency performed in the equatorial region on L ∼ 4 contain emissions of which a typical (f − t) spectrogram shows a pattern of stripes related to high-order proton cyclotron harmonics. These emissions have previously been attributed to the instability of ring-type distribution of energetic protons, based on the analysis of frequency dependence of the growth rate. In the present study we underline that the growth rate varies both in sign and in magnitude along the wave packet path; thus, the wave amplification, but not the growth rate, is the crucial quantity determining the observed spectrum. We develop a consistent model of the observed wave phenomenon using a smooth distribution function of energetic protons, which depends on particle energy W and equatorial pitch angle and, for a given W, tends to zero as → 0 and has maximum at = ∕2. The growth rate calculated with this distribution function differs essentially from that previously reported; in particular, has a negative minimum, instead of a positive maximum, at exact cyclotron resonances, that is, when the wave frequency is equal to a multiple of proton cyclotron frequency: = nΩ. This property is a key to understanding the peculiarity of the observed spectrum which consists in that often, although not always, the spectral intensity below the exact cyclotron harmonic is much higher than above. We have calculated the net amplification for a 3-D set of the wave packets. Assuming that the process of wave excitation is stationary and applying appropriate boundary conditions, we show that our model reproduces the wave phenomenon in outline.