Minimum variance projection for direct measurements of power-law spectra in the wavenumber domain

Minimum variance projection is widely used in geophysical and space plasma measurements to identify the wave propagation direction and the wavenumber of the wave fields. The advantage of the minimum variance projection is its ability to estimate the energy spectra directly in the wavenumber domain using only a limited number of spatial samplings. While the minimum variance projection is constructed for discrete signals in the data, we find that the minimum variance projection can reasonably reproduce the spectral slope of the power-law spectrum if the data represent continuous power-law signals. The spectral slope study using the minimum variance projection is tested against synthetic random data with a power-law spectrum. The method is applicable even for a small number of spatial samplings. Conversely, the spatial aliasing causes a flattening of the spectrum.


Introduction
Minimum variance estimator proposed by Capon (1969) (also referred to as the maximum likelihood method) offers a powerful tool for computing the wave energy directly in the wavevector domain using sensor array measurements and has been successfully applied to geophysical and space measurements.Examples are laboratory wave experiments for ocean wave breaking (Babanin et al., 2011), tidal interval waves (Dushaw and Worcester, 1998), oceanic Rossby waves (Zang and Wunsch, 1999), atmospheric layer imaging (Chen et al., 2014), seismic waves (Huang and Wu, 2006;Menon et al., 2014), and low-frequency waves in near-Earth space (magnetosheath) (Motschmann et al., 1996;Glassmeier et al., 2001).
The wave energy in Capon's method is evaluated from the covariance matrix of the measured data and the fluctuation model (referred to as the steering vector) as E mv (k) = 1 h † Rh .The covariance matrix is constructed from the measurement as R = ss t , where s = [s 1 , s 2 , . .., s n ] t denotes the state vector of the n-point measurement.The superscript t indicates the vector transposition.In wave studies, the steering vector is constructed as h = exp(ik • r 1 ), exp(ik • r 2 ), . .., exp(ik • r n ) using the scan wavevector k and the sensor coordinates r i (i = 1, 2, . .., n).The estimator is derived by minimizing the fluctuation variance projected onto the model geometry under a constraint of the unit gain for the signals (distortionless property), making use of the Lagrangian multiplier method.In other words, Capon's method is designed to detect the discrete signals in the measurement by enhancing the signal-tonoise ratio upon projection onto the model geometry.
Here we study the applicability of Capon's minimum variance projection in turbulence data analysis.Turbulence occurs in many geophysical flows and space plasmas, and the inertial range fluctuations are characterized by power-law Published by Copernicus Publications on behalf of the European Geosciences Union.
Y. Narita et al.: Minimum variance spectrum spectra.There are many signals in turbulence in the form of waves or eddies, and the number of sensors is by far smaller than that of signals.Our approach is numerical; we generate synthetic turbulence data, sample using a sensor array, and find the wave energy using Capon's method.We then compare Capon's minimum variance estimator spectra with that used in the synthetic data.
Taylor's frozen inflow hypothesis (Taylor, 1938) is not strictly valid even in the solar wind; thus, it is misleading to simply apply Taylor's hypothesis and to relabel the frequencies onto the stream-wise wavenumbers in interpreting the energy spectra from single-point measurements.The reason for this is the fact that the frequency broadening around the Doppler shift has substantial contributions, and the broadening can originate from not only the large-scale flow variations but also the existence of large-scale counterpropagating Alfvén waves (Narita, 2017).The spectral index can be the same between the frequency domain and the stream-wise wavenumber domain if the frequency slice of the spectrum at a given value of the wavenumber is Gaussian (Wilczek and Narita, 2012).In reality, however, there is no guarantee that the spectral index in the frequency domain is the same as that in the wavenumber domain since turbulent fluctuations can be non-Gaussian (coherent or intermittent) due to the growth or the decay of waves, the coherent structure formation, and the excitation of sideband waves.For this reason, there is a need to construct a method to determine the spectral indices directly in the wavenumber domain.

Numerical setup
The synthetic data are generated in the following steps.We set a mesh of one-dimensional spatial coordinates z n as z n = n z (n = 0, 1, . .., N − 1). (1) The number of mesh points is N = 2 12 = 4096 for the purpose of covering 3 orders of magnitude on the length scales or in the wavenumbers.The units of the spatial coordinate z are arbitrary.We set z = 1 km for a future application to in situ measurements in near-Earth space.The wavenumber mesh is constructed from the spatial coordinate mesh using the grid distance k = 2π z N −1 .
A power-law spectrum is defined at the mesh points of the wavenumbers as The spectral index α is a free parameter and is chosen in the range of typical values of turbulence energy spectra: The coefficient C is normalized to the variance of the fluctuations σ 2 through the relation where or We set the mean field as constant; B 0 = 1 nT and the standard deviation a σ = 0.1 nT (10 % fluctuation level).
The spectral amplitudes are computed in the wavenumber domain The wave phase φ(k) is chosen as uniformly and randomly distributed in the range 0 < φ ≤ 2π .The reality condition is used in the second half of the spectral data by taking the complex conjugate (denoted by the dagger † ): The fluctuations are obtained in the spatial coordinates using the inverse fast Fourier transform (FFT): The fluctuations in the spatial coordinates are displayed in the left panels of Fig. 1 for spectral indices of α = {−1.0,−1.2, . .., −2.0} (from top to bottom).
To double check the computation, the forward Fourier transform is applied to the spatial fluctuation data, and the energy spectra, which agree with the generating spectra, are obtained in the wavenumber domain (right panels in Fig. 1).The fluctuation data are defined in the spatial coordinates and do not contain temporal evolution.

Minimum variance spectra
The synthetic data are sampled at the positions of the sensor array, and the energy spectra are evaluated from the sampled data as follows.An array of four sensors is set with a sensor distance of = z = 1 km.The sensors are located on the linear mesh points.The fluctuation data are retrieved at the four sensor positions, and the energy spectrum is computed using the minimum variance estimator: The mean field B 0 = 1.0 nT is subtracted upon the minimum variance estimator.Isotropic level = 10 −8 nT 2 is added to the covariance matrix R to guarantee the matrix inversion.
The added isotropic level is well below the nominal sensor noise level of the fluxgate magnetometers.The energy spectra are averaged over 1000 realizations of the field measurements by shifting the sensor array in the positive z direction by a size of the sensor array, spanning from the first sensor to the last one (4 km).
The dimension of the averaged energy spectrum is converted from square amplitudes (nT 2 ) into an energy density in the spectral domain (nT 2 km rad −1 ) by dividing the square amplitude by the wavenumber interval for the sensor size, k = 2π L (where L = 4 ).The minimum variance estimator is designed to return a correct value of the spectral energy for a discrete signal in incoherent data.For multiple signals like turbulence, the spectral energy obtained by the minimum variance estimator is substantially diminished, and a correction factor is needed to normalize such that the integration of the spectrum over the wavenumbers is the same as the variance of the fluctuations.
The minimum variance estimator is tested under different spatial samplings: n = {2, 3, 4, 5, 10, 20, 50, 100} for a spectral slope of α = 2.The results are displayed in Fig. 2. The spectral energies are shifted by a factor of 10 −2 from one to another.The estimated spectral curves are smooth and deviated from the power law for a smaller number of spatial samplings (samplings of five points or fewer).The spectral curves exhibit fine structures and a well-established power law extending to lower and higher wavenumbers.

Spectral slopes
The major discovery from the numerical test for the minimum variance estimator is that the power-law indices or spectral slopes can reasonably be reproduced in a limited range of wavenumbers if the true spectrum follows a power law.Moreover, this ability to detect the true spectral slope is not lost, even using only a few spatial sampling points.
For comparison, we present two case studies for 4-point measurements (Fig. 3) and 20-point measurements (Fig. 4).The true spectra follow a power law with an index varying from −1 to −4.There is no spatial aliasing (see below), i.e., there is no signal with a wavelength below the sensor separation distance.The spectra are estimated using the minimum variance projection up to the Nyquist wavenumber k ny 3.1 rad km −1 .The true spectral slopes are reproduced at wavenumbers around 1 rad km −1 .
For the four-point measurements (Fig. 3), the lower limit of the useful wavenumber range for the power-law study depends on the true spectral curves: 0.8 rad km −1 for slopes −1 to −2, and 0.2 rad km −1 for a slope of −4.The upper limit of the useful wavenumber range is either up to the Nyquist wavenumber for slopes of −1 and −2 or about half of the Nyquist wavenumber for slopes of −3 and −4.
For the 20-point measurements (Fig. 4), the estimated spectra show well-established power-law curves and the spectral slopes agree with the true values.The range of useful wavenumbers extend to the Nyquist wavenumber in all the tested spectral slopes.The lower limit of the useful wavenumber range still depends on the true spectral slopes and varies from 0.8 rad km −1 for slopes of −1 and −2 to 0.2 rad km −1 for a slope of −4.

Aliasing problem
Since the spatial samplings are made at discrete points, the measurement of the wavenumber spectra is influenced by the spatial aliasing.That is, signals or fluctuations with a shorter wavelength than the sensor separation can affect the spectral measurements.In the case of the power-law spectra, the spatial aliasing affects the measurement through a flattening of the spectrum around the Nyquist wavenumber, k ny = π/d, where d represents a mean sensor separation distance in the array.Technically speaking, the minimum variance projection itself can compute the spectral energy beyond the Nyquist wavenumber limit, but the estimated spectra show either spikes due to the resonance for regularly spaced array sensors (we refer to it as regular aliasing) or a flat spectrum for irregularly spaced array sensor (we refer to it as irregular aliasing).
Examples of the regular and irregular aliasing-affected spectra are displayed in Figs. 5 and 6, respectively.
In Fig. 5, the case of regular aliasing, the signals are spatially down-sampled at the same distance of 10 km from one sensor to another, while the data are generated for the true power-law spectra up to a wavelength of 1 km with spectral slopes of −1, −2, −3, and −4.The Nyquist wavenumber is d 0.31 rad km −1 .Below the Nyquist wavenumber, the estimated spectra reproduce the true spectral slopes in a limited range of wavenumbers between 0.01 and 0.1 rad km −1 .Around the Nyquist wavenumber, the estimated spectra become flatter and increase toward the second-order harmonic of the Nyquist wavenumber at 2k ny .Beyond the Nyquist wavenumber, there is a series of spectral peaks at the higherorder harmonics of the Nyquist wavenumber, 2k ny , . .., 6k ny .The alias spikes at the harmonics of the Nyquist wavenum- ber can be seen both in the 4-sensor measurements (curves in black) and in the 20-sensor measurements (curves in gray).
In Fig. 6, the case of irregular aliasing, the signals are spatially down-sampled at a mean sensor separation distance of 10 km with a deviation of 30 % from the mean distance.The signals and the true spectra are the same as those used in the regular aliasing case.Again, the true spectral slopes can be reproduced in a relatively narrow range of the wavenumbers between 0.01 and 0.1 rad km −1 and the spectra show a flattening toward the Nyquist wavenumber.Beyond the Nyquist wavenumber, however, the alias spikes no longer appear and the spectral curves are nearly flat with a moderate fluctuation of the spectral energy.

Discussion and conclusions
It is encouraging that the minimum variance estimator is capable of not only identifying the wavenumbers in the discrete wave signals but also evaluating the spectral slope of the power-law spectrum in the continuous turbulence signals.The capability of the power-law estimation is still valid even using a small number of spatial samplings (typically below five spatial sampling points) at the cost of a limited wavenumber range and a smoothed shape of the spectral curve.Another finding from the presented numerical studies is that the spatial aliasing contributes in different ways.For regularly spaced (or equidistant) samplings, the minimum variance spectra show resonance peaks in the high wavenumber domain beyond the Nyquist wavenumber.For irregularly spaced samplings, the minimum variance spectra show no more resonance peak.
One does not need a frequency resolution since the technique is applied in the time domain.The strongest limitation is the fact that one has to assume that the length scales of the sensor array (both the separation distance and the total length) should be in the power-law spectrum range of the wavenumber domain.In the mean flow direction, one may use the time series data and Taylor's hypothesis to diagnose the existence of the power law in the corresponding wavenumber range (even though Taylor's hypothesis is not strictly valid).In the directions perpendicular to the mean flow, there is no means of diagnosing the existence of the power-law spectrum in the wavenumber domain.
The minimum variance estimator can be extended to a vectorial data set such as flow velocities or magnetic fields.In that case, the spectrum projected onto the wavevector domain is obtained as a matrix with the elements consisting of the covariance (represented in the wavevector domain) between different components of the measured vectorial field.
The minimum variance projection can be used for a variety of applications in space plasma observations.Examples are proposed below as a conclusion: one-dimensional application to magnetospheric plasma turbulence along the magnetic field line in the magnetotail using five THEMIS spacecraft (Angelopoulos, 2008); one-dimensional application to plasma turbulence in the inner heliosphere using multi-point measurements either along the magnetic field line or along the streamline by finding a suitable orbit configuration for Solar Orbiter (Müller et al., 2013), Solar Probe Plus (Fox et al., 2016), and BepiColombo cruise phase (Benkhoff et al., 2010); three-dimensional applications to local structures of plasma turbulence in space, e.g., shock-upstream and shock-downstream turbulence of the Earth magnetosphere, using tetrahedral spacecraft configuration of the Cluster mission (Escoubet et al., 2001) down to 100 km separation (corresponding to the spatial scale of ion gyroradius) and the Magnetospheric Multiscale (MMS) mission (Burch et al., 2016) down to 10 km separation (corresponding to electron gyroradius).
In individual cases, of course, the calibration of the minimum variance projection is needed as a pre-study because the spectral sensitivity (or gain) and the power-law range of wavenumbers might depend on the sensor array configuration.

Figure 1 .
Figure 1.Synthetic data (left panels) and the energy spectra used to generate the synthetic data (right panels).

Figure 2 .
Figure2.Minimum variance spectra for different numbers of spatial samplings for a spectral slope of −2.The spectra are shifted from one to another for the visualization purpose.

Figure 4 .
Figure 4. Minimum variance spectra for different slopes using 20 sampling points.Gray lines indicate the true spectral slope.

Figure 5 .
Figure 5. Minimum variance spectra with aliasing for regular spatial samplings (n = 4 in black and n = 20 in gray) with a distance of 10 km.

Figure 6 .
Figure6.Minimum variance spectra with aliasing for irregular spatial samplings (n = 4 in black and n = 20 in gray) with an irregular distance with a mean 10 km with 30 % deviation.
Narita et al.: Minimum variance spectrum Figure 3. Minimum variance spectra for different slopes using four sampling points.Gray lines indicate the true spectral slope.