Abstract

This paper presents a novel method for the automatic determination of the diagonal-loading level for robust adaptive beamforming on radar wind profilers. This method balances the degradation of the signal-to-interference ratio with that of the signal-to-noise ratio to maximize the detectability of the backscattered signals. Because radial wind velocities are usually estimated from the first moment of the spectrum of backscattered echoes, both the residual ground clutter and any increase in noise level degrade the detectability of atmospheric echoes. The proposed algorithm evaluates the power spectral density of the residual clutter and increased noise to determine the optimal diagonal-loading level by balancing these two factors. The results of numerical simulation show that, without the need to specify any user parameters, the proposed algorithm is stable and more effective at maximizing the signal-to-interference ratio than the conventional norm-constrained diagonal-loading approach. The stability and clutter suppression capability of the proposed algorithm are examined using data from the Program of the Antarctic Syowa Mesosphere–Stratosphere–Troposphere/Incoherent Scatter Radar.

1. Introduction

In radar wind profilers, suppression of ground clutter is critical in obtaining accurate measurements of wind velocities, especially when observing vertical air motion. Radial wind velocities are usually determined by taking the first moment of the spectrum of backscattered echoes, for example, using nonlinear spectral fitting methods (Sato and Woodman 1982). Vertical wind velocities are typically very small, so the spectral peak of turbulent echoes is close to zero Doppler velocity. Because ground clutter has a much stronger zero Doppler component in the spectrum, estimations of vertical wind velocities can contain serious errors.

Recently, phased arrays have become increasingly used as antenna structures for wind profilers. Phased array antennas have a number of advantages compared with mechanically scanned antennas, for example, rapid beam steering. In particular, by dividing the array into multiple subgroups, one can use adaptive beamforming techniques to mitigate clutter contamination. The ground clutter in wind profiler data is usually removed through spectral filtering, for example, by replacing the zero Doppler component by the average of both adjacent spectral bins. However, this procedure cannot distinguish between clutter and the actual atmospheric spectrum when it happens to be close to the zero Doppler component. Because adaptive beamforming uses the spatial distribution of targets, it is expected to improve the clutter rejection and hence the accuracy of vertical wind velocity measurements (Kamio and Sato 2004; Nishimura et al. 2012). The directionally constrained minimization of power (DCMP) algorithm (Takao et al. 1976), also known as the Capon beamformer (Capon 1969), is often used for adaptive clutter rejection because of its simplicity and effectiveness. The problem with using the DCMP algorithm, however, is that it can lead to uncontrollable changes in the main beam pattern. Thus, the reception beam direction might be changed considerably to suppress clutter, leading to severe degradation in the signal-to-noise ratio (SNR) and error in estimating the wind velocity using the Doppler beam swing method.

To control the amount of SNR degradation, the norm-constrained DCMP (NC-DCMP) algorithm was introduced. It is classified as a diagonal-loading beamformer with a diagonal-loading level determined by a designated norm constraint (Hudson 1981, 175–181). This method is known to be effective in actual observations from the middle and upper (MU) atmosphere radar in Shigaraki, Japan (Kamio and Sato 2004; Nishimura et al. 2012). More recently, Curtis et al. (2016) applied the same algorithm to observations from the National Weather Radar Testbed phased array radar. In the NC-DCMP algorithm, the norm constraint has clear physical meaning as the maximum permitted degradation in the SNR. Hence, it is usually set to a small value to maintain the noise floor level as low as possible and to enable the detection of weak atmospheric echoes. The selection of this value is empirical and not data adaptive; that is, the norm constraint must be determined in advance. However, if the norm constraint is too strict, then strong clutter might remain in the spectrum and cause severe bias to wind velocity estimations.

Recently, automatic determination of the diagonal-loading level has been extensively studied in pursuit of a more robust beamforming technique; a thorough review of these efforts was conducted by Du et al. (2009). They reported that previous algorithms have mainly focused on robustness against steering vector errors or a poorly estimated covariance matrix based on a small number of snapshots. For wind profilers, however, these problems are not critical. Indeed, the backscattered signals are assumed returned solely from the volume in the very sharp main lobe, and the direction of the center of this volume is determined by the transmitting beam pattern. Therefore, the steering vector for reception is the same as that for transmission. Similarly, the duration over which the covariance matrix is averaged can be relatively long (e.g., 1 min), yielding a sufficient number of snapshots for estimation of the covariance matrix with adequate accuracy. Instead, maintaining the white noise gain and beam directional errors as small as possible is critical because these attributes have greater impact on wind profilers, as they are expected to detect extremely weak signals. Currently, however, there are no suitable adaptive beamforming techniques that can automatically determine the optimal diagonal-loading level for such purposes.

In this paper, a novel method for automatic determination of the optimal diagonal-loading level for robust adaptive beamforming on radar wind profilers is presented. The algorithm developed here estimates the output powers of interference and noise simultaneously, and it automatically determines the optimal diagonal-loading level that balances the increase in the remaining clutter power with the increase in noise power to maximize the total detectability of the desired signal. Details of the proposed algorithm are described in section 3. The performance of this algorithm is evaluated and compared against both the standard DCMP algorithm and the conventional NC-DCMP algorithm in numerical experiments reported in section 4. The proposed algorithm is also applied to actual data from the Program of the Antarctic Syowa Mesosphere–Stratosphere–Troposphere/Incoherent Scatter Radar (PANSY; Sato et al. 2014) in section 5. Finally, its overall performance is summarized in section 6.

2. Basic methodologies

a. Nonadaptive beamforming

The narrowband beamforming output of a phased antenna array is generally written as

 
formula

where is a complex time series received by M spatially distributed receivers, is the weight vector, and Y is the synthesized output. The notation and represent the transposition and conjugate transposition of a matrix, respectively.

Nonadaptive beamforming uses the array manifold vector as its weight vector, where θ and ϕ are the zenith and azimuth angles, respectively, for an arbitrary direction in three-dimensional space. The ith component of the array manifold vector is written as

 
formula

where j is the imaginary unit, λ is the wavelength, is the location of the ith receiver, and is the radial unit vector in the given direction. Here, the azimuth angle is measured clockwise from north.

b. The DCMP algorithm

The DCMP algorithm is an adaptive beamforming method that minimizes the average output power under the constraint of the desired direction (Capon 1969; Takao et al. 1976). The DCMP algorithm can be written as a convex optimization problem,

 
formula

where is the covariance matrix and is the directional constraint. The solution to Eq. (3) is given as part of the description of the NC-DCMP algorithm.

Directional constraint is usually calculated using Eq. (2) with the desired direction , that is, , assuming uniform power directional gain for each channel. In applications of large atmospheric Doppler radars, however, the partially adaptive arrays should retain their main beam shape (Kamio and Sato 2004). In such configurations, a large main array is synthesized using nonadaptive beamforming to a single channel. A small number of antennas, called a subarray, are added to form a sidelobe canceller. Hence, the gain difference between the main and the subarray is usually sufficiently large to employ the following simple directional constraint: , where the first element corresponds to the main array (Curtis et al. 2016). When the abovementioned assumption fails, one can use the generalized gain-weighted directional constraint (Hashimoto et al. 2016). Each component of is weighted by the power directional gain of the corresponding antenna or subarray, which is written as

 
formula

where is the power directional gain of the ith element in the desired direction . Here, D contains the term to make , since is already normalized so that , as in Eq. (2). Note that when .

c. The NC-DCMP algorithm

The NC-DCMP algorithm is a modified DCMP that adds the following constraint to Eq. (3) (Hudson 1981, 175–181; Cox et al. 1987):

 
formula

where denotes the Euclidean norm and U is the norm constraint.

Although it is generally difficult to determine U, Kamio and Sato (2004) showed that U can be computed from the permissible SNR degradation caused by signal processing. As mentioned in the previous section, a partially adaptive array is preferable to a uniform-gain array for large atmospheric radars. Ideally, the desired signal is assumed unchanged by the DCMP algorithm in such nonuniform-gain configurations. As the output of the main array is synthesized by nonadaptive beamforming, it is dominant in determining the main beam shape. In addition, the contribution of the sidelobe canceller array to the main beam direction is kept small because of the gain weighting in the directional constraint. However, the average noise power is multiplied by because the primary noise source in the VHF band is the galactic noise that is random and independent of antenna gain. Therefore, the SNR degradation factor given by the standard DCMP algorithm is written as . Thus, to limit the SNR degradation to within , U can be set to . Term is usually set to a small value, such as (corresponding to ), which allows the algorithm to suppress the clutter at the cost of SNR degradation of less than 0.5 dB.

The optimal weight vector for the NC-DCMP algorithm can be obtained using the Lagrangian multiplier and diagonal-loading technique,

 
formula

where α is the diagonal-loading value and is the identity matrix. Note that gives the solution to the DCMP algorithm in Eq. (3). As decreases monotonically as α increases, the optimal α for the NC-DCMP algorithm is calculated as follows (Nishimura et al. 2012):

  1. Set α to a small value.

  2. Calculate using Eq. (6). If , then this is the solution.

  3. Otherwise, increase α and go to step 2.

The problem with the NC-DCMP algorithm is that the selection of the user-defined parameter is empirical. For example, Nishimura et al. (2012) used , 1.5, and 2.0, which correspond to the permissible SNR degradations having , 1.76, and 3 dB. As they pointed out, this parameter depends on the electromagnetic environment around the radar system: the background noise level, clutter characteristics, and power directionality of the antennas. In addition, these conditions also vary for each range, making selection of the norm constraint difficult. Therefore, an adaptive determination of the diagonal-loading value using observed data is desirable.

3. The proposed power balance algorithm

As described in the previous section, the NC-DCMP algorithm has a user-defined parameter to obtain the diagonal-loading level. This is not optimal in terms of the signal-to-interference ratio (SIR). Here, we propose a novel algorithm that automatically determines the optimal diagonal-loading level that balances the SNR and SIR degradations in dimensions of the power spectral density. In this section, we first show the derivations of the SIR and SNR degradations, followed by the formulation of the proposed algorithm.

a. Estimation of SNR and SIR degradations

Figure 1 shows the output power diagram of the desired signal , interference , and noise obtained using nonadaptive beamforming (), the standard DCMP algorithm (), and an intermediate diagonal-loading level . The total output power of the beamformer is written as follows [see also Eq. (3)]:

 
formula

where is calculated using Eq. (6) with a specific α. Note that is not a function of α. As mentioned in section 2c, the desired signal is considered unchanged by beamforming, especially when applied to a nonuniform-gain configuration. Hence, we assume

 
formula
Fig. 1.

Output power diagram for nonadaptive beamforming, the standard DCMP algorithm, and an intermediate diagonal-loading value α. Hatched area indicates the SIR degradation compared with the standard DCMP algorithm, and dotted area denotes the SNR degradation compared with nonadaptive beamforming.

Fig. 1.

Output power diagram for nonadaptive beamforming, the standard DCMP algorithm, and an intermediate diagonal-loading value α. Hatched area indicates the SIR degradation compared with the standard DCMP algorithm, and dotted area denotes the SNR degradation compared with nonadaptive beamforming.

One way of maximizing the signal-to-interference-plus-noise ratio (SINR) is to minimize the total output power in Eq. (7) because the signal power is unchanged, as stated in Eq. (8). Obviously, this describes the DCMP algorithm. Although the DCMP algorithm maximizes the SINR, it is not always optimal in terms of signal detectability. As stated in section 1, wind profilers generally use the Doppler spectrum to estimate radial wind velocity. In such applications, it is clearly important to retain the detectability of the signal; the spectral peak of atmospheric echoes should be as high as possible above the noise floor level (Gage and Balsley 1978). According to Kamio and Sato (2004), the DCMP algorithm can cause a severe increase in the sidelobe level when the desired signal is strong. This markedly degrades the detectability of the desired signal because of high white noise gain. Thus, the standard DCMP algorithm is unsuitable for wind profiler applications.

To address this problem, the residual clutter power and increased noise power must be evaluated separately. In Fig. 1, the dotted area denotes the SNR degradation compared with nonadaptive beamforming, and the hatched area indicates the SIR degradation compared with the standard DCMP algorithm. Below, we explain the derivations of these two quantities.

1) SNR degradation

As shown in section 2c, the output noise power is proportional to the norm of the weight vector. Therefore, the noise power obtained using α can be written as

 
formula

Using Eq. (9), the SNR degradation compared with the nonadaptive beamforming can be written as

 
formula
 
formula

To calculate using Eq. (11), must be estimated in advance. This is equal to the average power spectral density of the noise with appropriate normalization. In this paper, the average power spectral density of the noise is estimated using the segment method (Petitdidier et al. 1997). The segment method first divides the spectrum into equal-sized segments and then calculates the average power of each segment. Subsequently, the noise power is estimated by selecting the lowest average power with a factor that corrects for underestimation. This factor is calculated using a statistical property of the noise fluctuation that follows a distribution.

2) SIR degradation

The SIR degradation compared with the DCMP algorithm can be written as

 
formula

where denotes the interference power obtained using α. In contrast to in Eq. (9), cannot be directly estimated. However, the difference in the total output power related to the DCMP algorithm can be used to estimate ,

 
formula

Note that Eq. (13) has been simplified using the assumption in Eq. (8) that the signals are the same for different values of α. Rewriting Eq. (13) yields

 
formula

The first two terms on the right-hand side of this equation correspond to the upper hatched area in Fig. 1, which is the difference in the total output power. The other two terms correspond to the lower hatched area in Fig. 1, which can be interpreted as the hidden SIR degradation canceled by the SNR improvement related to employing the larger α.

3) Effects of diagonal loading on the SNR and SIR degradations

Figure 2 shows the relationship between and for α from 10−3 to 101.5. The horizontal axis represents the SIR degradation compared with the DCMP algorithm , and the vertical axis denotes the SNR degradation compared with the nonadaptive beamforming . The color of the circles on the line represents the magnitude of the diagonal-loading value α. Note that this example is taken from actual observations of PANSY that is discussed further in section 5. As in Fig. 2, there is a trade-off between the SIR and SNR degradations. Since both have the same dimensions, the residual clutter power and the increased noise power degrade the SINR by equivalent amounts.

Fig. 2.

Example of the relationship between the SIR degradation factor (abscissa) and SNR degradation factor (ordinate) for various diagonal-loading values α. Gray circles indicate the magnitude of α.

Fig. 2.

Example of the relationship between the SIR degradation factor (abscissa) and SNR degradation factor (ordinate) for various diagonal-loading values α. Gray circles indicate the magnitude of α.

b. Formulation of the cost function

So far, we have reviewed the main problem of the DCMP algorithm for wind profiler applications and derived SNR and SIR degradations as functions of α. We now consider the optimal cost function for the power minimization problem that is suitable for wind profilers.To balance the clutter suppression capability against signal detectability, the proposed algorithm solves the following minimization problem:

 
formula

Here, it should be recalled that both and degrade signal detectability by equivalent amounts. Hence, the optimal solution for Eq. (15) is the point on the curve in Fig. 2 that minimizes the distance to the origin (labeled as “optimal”). Using this optimal diagonal-loading level, the SNR and SIR degradations are balanced and the SINR in the spectral density ratio—that is, the SINDR—will be maximized. We call this method the power balance (PB) algorithm.

As shown in Eq. (15), we use the sum of squares of and . However, we can consider other cost functions, for example, the simple sum of these quantities; the cost function for which is written as . Again, this is equivalent to the standard DCMP algorithm, which can be naturally understood by rewriting using Eqs. (10) and (14) as follows:

 
formula
 
formula

where is a constant. Figure 3 shows examples of such cost function evaluations for and . The values used are the same as Fig. 2. The abscissa is the diagonal-loading value α and the ordinate is the cost function evaluation for the corresponding α. The solid and dashed lines denote the cost functions and , respectively. The color of the circles again indicates the magnitude of α. As shown in Fig. 3, the optimal solution for the standard DCMP algorithm is , because monotonically decreases as α decreases. In this case, the output SINR is maximized and clutter will be suppressed. However, increased noise power density might exceed the spectral peak of atmospheric echoes, and the signal could be completely lost. In contrast, the proposed PB algorithm evaluates both the amount of the suppressed clutter and the signal detectability separately. As a consequence, the cost function for the proposed algorithm reaches a minimum at around (labeled as “optimal”). At this point, signal detectability and clutter suppression are balanced and this eventually maximizes the SINDR.

Fig. 3.

Example of the relationship between various diagonal-loading values α and the cost functions of the PB algorithm, , and standard DCMP algorithm, .

Fig. 3.

Example of the relationship between various diagonal-loading values α and the cost functions of the PB algorithm, , and standard DCMP algorithm, .

c. Procedure for the power balance algorithm

The procedure for the proposed PB algorithm is as follows:

  1. Estimate the mean noise power using nonadaptive beamforming .

  2. Calculate the optimal weight of the standard DCMP algorithm , total output power , and noise power .

  3. Find the value of α that minimizes the cost function by solving Eq. (15).

  4. If α is greater than the minimum diagonal-loading value ε, then this is the optimal solution, otherwise use .

Here, we use a minimum diagonal-loading value ε because the loss of signal power caused by the DCMP algorithm cannot be ignored in some situations, for example, when the desired signal is very strong and clutter is absent. In these cases, the assumption about the signal power in Eq. (8) will be violated, resulting in an incorrect estimation of the SIR degradation in Eq. (14). To reduce the loss of signal power to a negligible amount, ε can be calculated as follows:

 
formula

which is equivalent to the SNR degradation given by the DCMP algorithm. Supposing there is only the desired signal and noise, the loss of signal power will be at least larger than ε because the total output power of the standard DCMP algorithm must always be smaller than the nonadaptive beamforming. For this reason, we define the minimum diagonal-loading value of to compensate for this loss.

Now, we consider the computational cost of the proposed PB algorithm. The method is in the class of nonlinear least squares problems. Hence, Eq. (15) can be efficiently solved using the Levenberg–Marquardt method (Marquardt 1963). Because Eq. (6) is evaluated at every step to find the optimal α, the difference in the computational cost is proportional to the number of evaluations of Eq. (6). In this paper, the inverse of the diagonally loaded covariance matrix is calculated using eigendecomposition (Hudson 1981, 175–181),

 
formula

where is the eigendecomposition of the covariance matrix, is the ith eigenvalue, and is the ith eigenvector. Note that in the VHF band, galactic noise power is usually sufficiently large to make positive definite, that is, all . Equation (19) enables rapid computation of Eq. (6) because eigendecomposition is required only once for each covariance matrix. Thus, the increase in computational complexity for the proposed PB algorithm is small.

Here, we provide an example based on actual observations, which is described in more detail in section 5. In this example, each signal block has eight channels, 1024 time samples, and 157 range samples. Because the sampling interval is 51.2 ms, each block must be processed within about 52 s. In our implementation, the proposed PB algorithm required 2–3 times as many cost function evaluations as the NC-DCMP algorithm. Consequently, the average computational durations were 5.40 s for the PB, 3.68 s for the NC-DCMP, and 2.74 s for the DCMP algorithms on a personal computer capable of performing 112 billion floating-point operations per second. All these computational times were sufficiently shorter than the data interval of 52 s. Therefore, we can conclude that the computational cost of the proposed PB algorithm is comparable to the conventional NC-DCMP algorithm and that it could be used in real-time applications.

In the following section, the performance of the proposed PB algorithm is compared with those of the standard DCMP algorithm and the conventional NC-DCMP algorithm using numerical simulations.

4. Numerical simulations

a. System model

In this simulation, the radar system is based on the MU radar (Hassenpflug et al. 2008). The MU radar consists of 475 three-element crossed-Yagi antennas. The antenna array is divided into 25 subarrays, each consisting of 19 antennas. In this simulation, only 23 subarrays are used, which is the same setting as for the observations used in section 5. The antenna arrangement of the MU radar is illustrated in Fig. 4. Other simulation parameters listed in Table 1 are the same as for the observations in section 5.

Fig. 4.

Antenna position and subarray configuration of the MU radar used for simulations and observations. White circles are antennas in main array, and black circles are for sidelobe canceller array.

Fig. 4.

Antenna position and subarray configuration of the MU radar used for simulations and observations. White circles are antennas in main array, and black circles are for sidelobe canceller array.

Table 1.

Parameters for simulations and observations made on 2 Jul 2015 by the MU radar.

Parameters for simulations and observations made on 2 Jul 2015 by the MU radar.
Parameters for simulations and observations made on 2 Jul 2015 by the MU radar.

As mentioned in section 2b, a nonuniform-gain configuration is preferable to retain the shape of the main lobe. To confirm this is the case, both uniform-gain and nonuniform-gain configurations are considered in this simulation. In both configurations, the entire array is divided into six groups. For the uniform-gain configuration, each of the nearest three groups arranged in a regular triangle are synthesized in phase.

For the nonuniform-gain configuration, the output from all hexagons is synthesized by nonadaptive beamforming. Additionally, five subarrays added as sidelobe cancellers are indicated by black circles in Fig. 4. Figure 5 shows the power directionality of the main array and one of the elements of the subarray in the azimuth section at .

Fig. 5.

Power directionality pattern of the main array (Main) and sidelobe canceller (SC) array of the MU radar used for simulations and observations. This is the section at an azimuth angle of .

Fig. 5.

Power directionality pattern of the main array (Main) and sidelobe canceller (SC) array of the MU radar used for simulations and observations. This is the section at an azimuth angle of .

b. Signal generation

In this simulation, there are three types of signal: atmospheric echoes, noise, and stationary clutter. The detailed procedures for generating these signals are described below.

1) Atmospheric echoes

Atmospheric echoes are generated using the atmospheric backscatter simulator developed by Holdsworth and Reid (1995). Here, we present the simulation parameters used for generating atmospheric echoes and a brief explanation of this methodology. At the beginning of a simulation, a finite number of scatterers are randomly placed in the target radar volume, which is defined by a cylinder enclosing the effective beamwidth and the target range gate. The time series measured at each receiver becomes the sum of the signals from all scatterers. The amplitude of the signal from each scatterer is determined from the random reflectivity ratio assigned to each scatterer as well as the weighting function for the position of the scatterer defined by the beam pattern and the range weighting function. The phase is determined by the distance between each receiver and each scatterer. For every sample time, the position of each scatterer is updated together with the background wind and the turbulent motion generated by mutual interaction among the scatterers. Turbulent motion is calculated using the random turbulent motion vector assigned to each scatterer. Scatterers that have moved outside of the enclosing volume are reentered from the opposite side with new values for random reflectivity and their turbulent motion vector. Table 2 lists the simulation parameters used to generate atmospheric echoes. Note that the signal-to-noise density ratio (SNDR) is defined as the spectral peak of atmospheric echoes compared with the noise floor level. The time series generated using this simulator for each receiver is written as , where k is the sampling time index.

Table 2.

Simulation parameters used for generating atmospheric backscatter signals.

Simulation parameters used for generating atmospheric backscatter signals.
Simulation parameters used for generating atmospheric backscatter signals.

2) Noise

Noise is modeled as a complex random number having real and imaginary parts that follow a zero-mean, unit-variance normal distribution. The resultant noise sequence is normalized by to make the average noise power .

3) Ground clutter

Ground clutter is modeled by point targets at low elevation angles with random directions. In this simulation, five point sources are generated at directions determined by a uniformly distributed random number in the range for the zenith angle and for the azimuth angle. Here, we use five sources because this is the maximum number of clutter sources that can be handled within the system in this simulation. The complex time series of the received signal at each receiver can be modeled as

 
formula

where is the total power from all clutter sources; is the direction of the jth clutter source; and and are the array manifold vector and the directionality gain of the ith receiver, respectively. Note that the right-hand side of Eq. (20) is a complex constant that depends only on the incident angles , because each ground clutter signal is modeled as a stationary source and its location is independent of time. The clutter power is selected such that the total signal-to-interference density ratio (SIDR) obtained by nonadaptive beamforming is equal to the designated value. Here, the SIDR is defined as the difference between the atmospheric and clutter spectral peaks. The periodogram of the clutter signal with eight-time incoherent integration is calculated as

 
formula

where ω is the Doppler frequency, denotes the discrete-time Fourier transform, is the optimal weight vector calculated from the covariance matrix using a diagonal-loading value α, is the sample index in the mth periodogram, and is the length of each periodogram. Then, is determined such that the interference-to-noise density ratio (INDR) satisfies the following equation:

 
formula

where is the zero Doppler frequency component of the periodogram . Here, INDR is defined as the peak of the clutter spectrum compared with the noise floor level.

c. Signal processing

For both uniform- and nonuniform-gain configurations, the DCMP algorithm, NC-DCMP algorithm, and proposed PB algorithm are applied to the received signal , defined by

 
formula

The sample covariance matrix is calculated using snapshots around the sampling time index k,

 
formula

The sampling interval is ms, which is equivalent to time averaging over approximately 26 s. The permissible SNR degradation of the NC-DCMP algorithm is set to 0.5 dB, which corresponds to the norm constraint . As listed in Table 2, the beam direction is . All antennas have the same element gain function; the gain weighting coefficients of the nonuniform-gain configuration are proportional to the number of antennas in the main and subarrays, that is, 874 and 19, respectively. Therefore, we use and .

d. Performance evaluation method

Once the optimal weight vector has been calculated, the performance of each method for both configurations is measured by the SINDR, beam directional error, and quality index. To obtain statistics for both the uniform- and nonuniform-gain configurations, 100 Monte Carlo simulations were executed for two SNDR cases: 20 and 40 dB. In each simulation, the input SIDR was increased from to 0 dB at intervals of 10 dB.

The SINDR is the SINR for the spectral density, which is defined as

 
formula

where is the Doppler frequency corresponding to the peak of the synthesized periodogram of the atmospheric echo obtained by α and is the noise power with α. Note that is the averaged periodogram calculated in the same manner as the clutter signal using Eq. (21).

To calculate the beam directional error, the beam direction is first calculated by determining the direction of the most significant peak of the power directionality pattern with the optimal weight vector,

 
formula

The beam directional error is then calculated by taking the average of the angles between the two vectors and ,

 
formula

where denotes the inner product of the vectors and .

The quality index Q is defined as the product of the SNDR and SIDR degradation indices L and Z,

 
formula

Larger values () mean better performance. The SNDR degradation index L is defined as

 
formula

where the numerator and denominator correspond to the SNDR of the nonadaptive beamforming and a particular α, respectively. The value of L is always less than 1 and larger values of L mean less SNR degradation. The SIDR degradation index Z is defined as

 
formula

which is the ratio between the SIDR at a particular α and the SNDR of the nonadaptive beamforming. A value of means that the clutter is sufficiently suppressed and that its peak is below the noise floor level. Values of do not improve the detectability of the signal; therefore, it is clipped at 1.

e. Results and discussion of simulation

Figures 6a–f show the average SINDRs and beam directional errors for the DCMP algorithm, conventional NC-DCMP algorithm, and proposed PB algorithm for both uniform- and nonuniform-gain configurations. The input SIDRs range from to 0 dB, while the input SNDR is 5 dB for Figs. 6a and 6b, 20 dB for Figs. 6c and 6d, and 40 dB for Figs. 6e and 6f. Figure 7 shows examples of the reception beam patterns for the DCMP, NC-DCMP, and PB algorithms when the input SNDR and SIDR are 40 and 0 dB, respectively. The azimuth angle for Fig. 7 is .

Fig. 6.

(a),(c),(e) Average SINDRs for the input SNDR of 5, 20, and 40 dB. (b),(d),(f) Average beam directional error compared with the nonadaptive beamforming for the input SNDR of 5, 20. and 40 dB in the simulation. Dotted lines are for the DCMP algorithm (DCMP), dashed lines for the NC-DCMP algorithm (NC), and solid lines for the proposed PB algorithm (PB). Lines without markers are for the uniform-gain configuration (U), and lines with markers are for the nonuniform-gain configurations (N).

Fig. 6.

(a),(c),(e) Average SINDRs for the input SNDR of 5, 20, and 40 dB. (b),(d),(f) Average beam directional error compared with the nonadaptive beamforming for the input SNDR of 5, 20. and 40 dB in the simulation. Dotted lines are for the DCMP algorithm (DCMP), dashed lines for the NC-DCMP algorithm (NC), and solid lines for the proposed PB algorithm (PB). Lines without markers are for the uniform-gain configuration (U), and lines with markers are for the nonuniform-gain configurations (N).

Fig. 7.

Example of the reception beam pattern synthesized by the DCMP algorithm (DCMP), NC-DCMP algorithm (NC), and proposed PB algorithm (PB) when the input SNDR is 40 dB. This is the section at an azimuth angle of .

Fig. 7.

Example of the reception beam pattern synthesized by the DCMP algorithm (DCMP), NC-DCMP algorithm (NC), and proposed PB algorithm (PB) when the input SNDR is 40 dB. This is the section at an azimuth angle of .

As shown by Figs. 6b, 6d, and 6f, the nonuniform-gain configuration mostly gives a better result than the uniform-gain configuration. In particular, Fig. 6d shows that the beam directional errors are at least about for the NC-DCMP algorithm with the uniform-gain configuration, while less than for the nonuniform-gain configuration. In vertical wind measurements, beam directional errors can cause leakage of the horizontal wind speed into vertical wind components, which can cause severe estimation error. These correspond to vertical wind speed errors of about 0.21 and 0.06 m s−1, respectively, when the horizontal wind speed is 40 m s−1, which is a typical value of the upper-tropospheric jet. Thus, the improvement of wind estimation accuracy by employing the nonuniform-gain configuration is about 0.15 m s−1. This is considerable improvement because the vertical wind speed is usually very small, for example, mostly less than about 0.4 m s−1 in observations from PANSY (Sato et al. 2014).

The beam directional error has a serious effect, especially on the estimation of the mean vertical wind, which is of the order of 0.01 m s−1 after the averaging of about half an hour. It is an important physical parameter in the discussion of atmospheric dynamics, such as convection, gravity waves, and ageostrophic motions associated with weather disturbances. To estimate it, an error of at most 0.1 m s−1 is desirable for each spectrum before averaging (e.g., Fukao et al. 1991; Nishimura et al. 2012).

Thus, the nonuniform-gain configurations are desirable for the Doppler beam swing method, although these errors could be corrected by calculating the actual beam direction using each weight vector with Eq. (26). These results support our discussion of the advantages of a nonuniform-gain configuration in section 2c; the noise power increase and loss of signal power are kept small through the gain weighting. As shown in Fig. 6f, the beam directional errors become more severe for the uniform-gain configurations when the SNDR is 40 dB, which causes greater inaccuracy in the estimations of wind velocity. Thus, we discuss the results for only the nonuniform-gain configuration hereafter.

As shown in Fig. 6c, the DCMP algorithm performs similarly to the PB algorithm, although its SINDR is about 0.8 dB lower for SIDRs below dB. For SIDRs above dB, the PB algorithm has lower SINDRs than the DCMP algorithm, namely, about 1 dB. This is because clutter is not sufficiently suppressed in this region, even by the DCMP algorithm. As shown by Kamio and Sato (2004), the depth of a null formed by the DCMP algorithm depends on the strength of clutter; that is, it is difficult to suppress weak clutter. As the PB algorithm measures SIR degradation using the difference calculated with the DCMP algorithm, the residual clutter power estimated by the PB algorithm tends to be larger than the DCMP algorithm. Therefore, if clutter is left by the DCMP algorithm, the PB algorithm also fails to estimate the actual clutter power, resulting in higher residual clutter power than for the DCMP algorithm. In contrast, the output SINDRs of the NC-DCMP algorithm are inferior to the PB and DCMP algorithms because of the predefined norm constraint; that is, is too strict in this case. As a result, the NC-DCMP algorithm shows the smallest beam directional errors, as illustrated in Fig. 6d. Although the output SINDRs would be higher if a larger norm constraint were set, no relevant criterion exists for determining the optimal norm constraint for the NC-DCMP algorithm, as mentioned in section 2c.

Although the SINDR of the DCMP algorithm is comparable to the PB algorithm in the previous case, it degrades markedly when the input SNDR is 40 dB, as shown in Fig. 6e. The SINDR of the DCMP algorithm is on average about 10 dB lower than that of the PB algorithm. This is because the loss of signal power is not negligible for the DCMP algorithm; that is, the shape of the beam pattern is significantly changed, as shown in Fig. 7. The NC-DCMP algorithm also shows lower SINDRs than the PB algorithm, except for an SIDR of 0 dB. Here, the NC-DCMP algorithm shows the best SINDR because the current norm constraint is closer to the optimum than in the previous case. Even in such situations, however, the PB algorithm uses the minimum diagonal-loading value to preserve the main beam shape, as described in section 3c, and the SNDR degradations are limited compared with the DCMP algorithm.

On the other hand, as shown in Fig. 6a, the output SINDRs of the DCMP algorithm are the highest for most of the input SIDRs when the input SNDR is 5 dB. This is because when the desired signal power is very low, the loss of the desired signal also becomes very small, even when using the DCMP algorithm. However, since the input SNDR cannot be known in advance, it is difficult to determine whether it is safe to use the DCMP algorithm.

These results are summarized in Fig. 8 for the nonuniform-gain configuration showing, from left to right, the SIDR degradation index Z, SNDR degradation index L, and performance index Q for the PB (PB), NC-DCMP (NC), and DCMP algorithms. The central line of each box shows the median; marks are the means; upper and lower edges of the box are the first and third quartiles, respectively; and whiskers are the upper and lower interquartile ranges multiplied by 1.5. Cases from all SIDR variations and all SNDRs are averaged.

Fig. 8.

Clutter suppression capability Z, SNR loss L, and quality index Q for the DCMP algorithm (DCMP), NC-DCMP algorithm (NC), and proposed PB algorithm (PB) in the simulation for nonuniform-gain configuration.

Fig. 8.

Clutter suppression capability Z, SNR loss L, and quality index Q for the DCMP algorithm (DCMP), NC-DCMP algorithm (NC), and proposed PB algorithm (PB) in the simulation for nonuniform-gain configuration.

As shown in the left panel of Fig. 8, the SIDR degradation of the PB algorithm is on average the smallest of all the methods. Although the clutter suppression capability of the DCMP algorithm is considered the highest, Z is not the highest because the signal degradations are dominant, as mentioned above. On the other hand, as shown in the middle panel, the SNDR degradation of the NC-DCMP algorithm is the smallest of all the methods. The SNDR degradation of the PB algorithm is on average about 3.9 dB larger to improve the clutter suppression capability by about 17.4 dB. As a result, the PB algorithm shows the best performance on average, as shown in the right panel of Fig. 8.

From the abovementioned discussion, we conclude that the proposed PB algorithm has desirable characteristics for wind profilers and that it could replace the conventional NC-DCMP algorithm. The proposed PB algorithm shows the best compromise compared with the DCMP and NC-DCMP algorithms—that is, better robustness than the DCMP algorithm—as the SNDR varies and higher SINDR than the conventional NC-DCMP algorithm as the SIDR varies, when applied to nonuniform-gain configurations. In addition, its diagonal-loading level is automatically calculated to balance the SNR and SIR degradations, making this algorithm free of the need for user-defined parameters.

5. Application to radar observations

a. Observations

Two observation datasets are used in this section. The first was obtained on 2 July 2015 using the MU radar in Shigaraki, Japan (Hassenpflug et al. 2008). The settings are the same as those used in the simulation described in section 4, except that only the nonuniform-gain configuration is considered because it shows much smaller beam directional errors. Here, the north beam directed to is used. Observations from 0900 to 1000 UTC are used for averaging.

The second dataset was obtained on 20 March 2015 using PANSY at Syowa station, Antarctica (Sato et al. 2014). PANSY consists of 1045 three-element crossed-Yagi antennas. The antenna array is divided into 55 subarrays, each of which consists of 19 antennas. In this observation dataset, two subarrays were not used, reducing the total number of subarrays to 53. To form the nonuniform-gain adaptive array, seven subarrays were used as sidelobe cancellers, which are indicated by hexagonal frames in Fig. 9. The other 46 subarrays, indicated by white circles, were synthesized by nonadaptive beamforming. Other observation parameters are shown in Table 3. Again, the beam direction of is used. Observations from 1236 to 1315 UTC are used for averaging.

Fig. 9.

Antenna position and subarray configuration of PANSY used in the observations. Each set of antennas surrounded by a hexagonal frame indicates a single channel of the sidelobe canceller subarray. All other antennas with white circles are synthesized in phase.

Fig. 9.

Antenna position and subarray configuration of PANSY used in the observations. Each set of antennas surrounded by a hexagonal frame indicates a single channel of the sidelobe canceller subarray. All other antennas with white circles are synthesized in phase.

Table 3.

Parameters for observations made on 20 March 2015 by PANSY.

Parameters for observations made on 20 March 2015 by PANSY.
Parameters for observations made on 20 March 2015 by PANSY.

b. Signal processing

The standard DCMP algorithm, NC-DCMP algorithm, and proposed PB algorithm are applied to a nonuniform-gain configuration with eight receivers, as described in section 4a. The permissible SNR loss for the NC-DCMP algorithm is set to 0.5 dB.

The performance of the beamformer is evaluated using a periodogram with frequency bins. Before incoherent integration, the noise floor level is corrected using the average squared norm of the weight vector in a range-by-range manner. Each periodogram is written as

 
formula

where is the sampling time index. As discussed in section 2c, the increase in the noise floor level in the periodogram is proportional to the squared norm of the weight vector. Hence, is divided by the squared norm of the weight vector to obtain the normalized periodogram ,

 
formula

Note that this correction is not only required to compare the results of different adaptive beamforming algorithms but also to ensure across range continuity of the noise power; otherwise, the noise floor level might have a different bias at each range location.

The performance indices Z, L, and Q are calculated in the same manner as in section 4d.

c. Results and discussion of application to radar observations

Figure 10 shows the indices Z, L, and Q for the PB (PB), NC-DCMP (NC), and DCMP algorithms for the first dataset from the MU radar. The styles are the same as in Fig. 8. We used 5472 periodograms to calculate the statistical results.

Fig. 10.

Clutter suppression capability Z, SNR loss L, and quality index Q for the DCMP algorithm (DCMP), NC-DCMP algorithm (NC), and proposed PB algorithm (PB) in the observation of the MU radar.

Fig. 10.

Clutter suppression capability Z, SNR loss L, and quality index Q for the DCMP algorithm (DCMP), NC-DCMP algorithm (NC), and proposed PB algorithm (PB) in the observation of the MU radar.

As shown in Fig. 10, the result from the observation is consistent with that of the simulation in Fig. 8. The PB algorithm indicates the highest average clutter suppression capability at the cost of limited SNDR degradation; that is, it is about 3 dB larger than the NC-DCMP algorithm. The resultant performance index Q of the PB algorithm is, on average, the highest among all three methods. Average values of Q in Fig. 10 are , , and for the PB, NC-DCMP, and DCMP algorithms, respectively.

Figures 11a and 11b show the Doppler spectra at 6.0 and 4.2 km, respectively, from the second dataset from PANSY. These are after 60-time incoherent integration, which corresponds to about 40 min. The horizontal axis is the Doppler velocity and the vertical axis is the power spectral density. The thin solid line denotes the nonadaptive beamforming (NA), the dotted line denotes the standard DCMP algorithm (DCMP), the dashed line denotes the NC-DCMP algorithm (NC), and the thick solid line with marks denotes the proposed PB algorithm (PB).

Fig. 11.

Doppler spectra averaged for about 40 min at (a) 6.0 and (b) 4.2 km obtained by PANSY using the nonadaptive beamforming (NA), DCMP algorithm (DCMP), NC-DCMP algorithm (NC), and proposed PB algorithm (PB).

Fig. 11.

Doppler spectra averaged for about 40 min at (a) 6.0 and (b) 4.2 km obtained by PANSY using the nonadaptive beamforming (NA), DCMP algorithm (DCMP), NC-DCMP algorithm (NC), and proposed PB algorithm (PB).

As seen in Fig. 11a, strong ground clutter is left in the periodogram obtained from the NC-DCMP algorithm at 6.0 km. However, this ground clutter is sufficiently suppressed by the PB and DCMP algorithms at a cost of about 1-dB additional SNDR degradation compared with the NC-DCMP algorithm. The average α for this range is for the NC-DCMP algorithm and for the PB algorithm. Such a large difference in α values implies that the ground clutter in this range is difficult to suppress using the small norm constraint of .

In contrast to the 6.0-km case, Fig. 11b indicates that the ground clutter is suppressed by all three methods at 4.2 km. This is probably because the directionality pattern response of the subarray at the incident angle of the ground clutter is higher than in the previous case. The diagonal-loading values support this presumption because they are much smaller, that is, and .

Now we compare the results for the PB and DCMP algorithms. The SNDR degradation for the PB algorithm is about 0.4 dB less than the DCMP algorithm, as shown in the zoomed portions of Figs. 11a and 11b.

On the other hand, the clutter suppression capabilities of these two methods are the same in both figures. This agrees with our simulation results in section 4; for example, when the SNDR is about 20 dB and the SIDR is around dB, the PB and DCMP algorithms work similarly. However, we have also shown in our simulation that the DCMP algorithm can cause severe SNR degradation with high sidelobes when the input SNDR is high. Therefore, the PB algorithm is considered a better solution, even when SNR degradation caused by DCMP algorithm is not severe, because the PB algorithm can prevent an unpredictable noise power increase.

6. Summary and conclusions

This paper describes a novel method for the automatic determination of the diagonal-loading level for robust adaptive beamforming on radar wind profilers by balancing the SIDR and SNDR degradations to maximize the detectability of the desired signals. The proposed PB algorithm evaluates the residual clutter power and increased noise power in the power spectral density, making the algorithm suitable for applications dealing with extremely weak signals. The proposed method also shows robustness against high SNDRs when the performance of the standard DCMP algorithm deteriorates with high sidelobes. The algorithm includes a nonlinear least squares problem that increases its complexity compared with the conventional NC-DCMP algorithm. However, the computational complexity remains sufficiently small to be applied in real-time applications.

In section 4, the performance of the proposed algorithm was examined using numerical simulations. Without the need to specify any user-defined parameters, the proposed algorithm achieved, on average, both better clutter suppression and better robustness against the input SNDR variation compared with the conventional NC-DCMP and DCMP algorithms. In section 5, the proposed algorithm was applied to observations from the MU and PANSY radars. The proposed algorithm sufficiently suppressed ground clutter with a smaller noise floor increase than the standard DCMP algorithm, even in ranges where the conventional NC-DCMP algorithm failed. In addition, the results of the statistics are consistent with the simulation, which supports the assertion made in section 4 that the proposed PB algorithm is the best compromise against the SIDR and SNDR variations in real applications. From these results, we conclude that the proposed PB algorithm is suitable for wind profilers and that it could replace the conventional NC-DCMP algorithm.

Acknowledgments

Shigaraki MU Observatory is operated by the Research Institute for Sustainable Humanosphere, Kyoto University. PANSY is a multi-institutional project with a core team at both the University of Tokyo and the National Institute of Polar Research. PANSY is operated by the Japanese Antarctic Research Expedition.

REFERENCES

REFERENCES
Capon
,
J.
,
1969
:
High-resolution frequency-wavenumber spectrum analysis
.
Proc. IEEE
,
57
,
1408
1418
, doi:.
Cox
,
H.
,
R. M.
Zeskind
, and
M. M.
Owen
,
1987
:
Robust adaptive beamforming
.
IEEE Trans. Acoust. Speech Signal Process.
,
35
,
1365
1376
, doi:.
Curtis
,
C. D.
,
M.
Yeary
, and
J. L.
Lake
,
2016
:
Adaptive nullforming to mitigate ground clutter on the National Weather Radar Testbed Phased Array Radar
.
IEEE Trans. Geosci. Remote Sens.
,
54
,
1282
1291
, doi:.
Du
,
L.
,
T.
Yardibi
,
J.
Li
, and
P.
Stoica
,
2009
:
Review of user parameter-free robust adaptive beamforming algorithms
.
Digital Signal Process.
,
19
,
567
582
, doi:.
Fukao
,
S.
,
M. F.
Larsen
,
M. D.
Yamanaka
,
H.
Furukawa
,
T.
Tsuda
, and
S.
Kato
,
1991
:
Observations of a reversal in long-term average vertical velocities near the jet stream wind maximum
.
Mon. Wea. Rev.
,
119
,
1479
1489
, doi:.
Gage
,
K. S.
, and
B. B.
Balsley
,
1978
:
Doppler radar probing of the clear atmosphere
.
Bull. Amer. Meteor. Soc.
,
59
,
1074
1093
, doi:.
Hashimoto
,
T.
,
K.
Nishimura
, and
T.
Sato
,
2016
:
Adaptive sidelobe cancellation technique for atmospheric radars containing arrays with nonuniform gain
.
IEICE Trans. Commun.
,
E99.B
,
2583
2591
, doi:.
Hassenpflug
,
G.
,
M.
Yamamoto
,
H.
Luce
, and
S.
Fukao
,
2008
:
Description and demonstration of the new Middle and Upper atmosphere Radar imaging system: 1-D, 2-D, and 3-D imaging of troposphere and stratosphere
.
Radio Sci.
,
43
,
RS2013
, doi:.
Holdsworth
,
D. A.
, and
I. M.
Reid
,
1995
:
A simple model of atmospheric radar backscatter: Description and application to the full correlation analysis of spaced antenna data
.
Radio Sci.
,
30
,
1263
1280
, doi:.
Hudson
,
J.
,
1981
: Adaptive Array Principles. IET Electromagnetic Waves Series, Vol. 11, IET, 268 pp., doi:.
Kamio
,
K.
, and
T.
Sato
,
2004
:
An adaptive sidelobe cancellation algorithm for high-gain antenna arrays
.
Electron. Commun. Japan
,
87
,
11
18
, doi:.
Marquardt
,
D. W.
,
1963
:
An algorithm for least-squares estimation of nonlinear parameters
.
J. Soc. Ind. Appl. Math.
,
11
,
431
441
, doi:.
Nishimura
,
K.
,
T.
Nakamura
,
T.
Sato
, and
K.
Sato
,
2012
:
Adaptive beamforming technique for accurate vertical wind measurements with multichannel MST radar
.
J. Atmos. Oceanic Technol.
,
29
,
1769
1775
, doi:.
Petitdidier
,
M.
,
A.
Sy
,
A.
Garrouste
, and
J.
Delcourt
,
1997
:
Statistical characteristics of the noise power spectral density in UHF and VHF wind profilers
.
Radio Sci.
,
32
,
1229
1247
, doi:.
Sato
,
K.
, and Coauthors
,
2014
:
Program of the Antarctic Syowa MST/IS radar (PANSY)
.
J. Atmos. Sol.-Terr. Phys.
,
118A
,
2
15
, doi:.
Sato
,
T.
, and
R.
Woodman
,
1982
:
Spectral parameter estimation of CAT radar echoes in the presence of fading clutter
.
Radio Sci.
,
17
,
817
826
, doi:.
Takao
,
K.
,
M.
Fujita
, and
T.
Nishi
,
1976
:
An adaptive antenna array under directional constraint
.
IEEE Trans. Antennas Propag.
,
24
,
662
669
, doi:.

Footnotes

© 2017 American Meteorological Society. For information regarding reuse of this content and general copyright information, consult the AMS Copyright Policy (www.ametsoc.org/PUBSReuseLicenses).