Compressed Symmetric Nested Arrays and Their Application for Direction-of-Arrival Estimation of Near-Field Sources

In this paper, a new sensor array geometry, called a compressed symmetric nested array (CSNA), is designed to increase the degrees of freedom in the near field. As its name suggests, a CSNA is constructed by getting rid of some elements from two identical nested arrays. The closed form expressions are also presented for the sensor locations and the largest degrees of freedom obtainable as a function of the total number of sensors. Furthermore, a novel DOA estimation method is proposed by utilizing the CSNA in the near field. By employing this new array geometry, our method can identify more sources than sensors. Compared with other existing methods, the proposed method achieves higher resolution because of increased array aperture. Simulation results are demonstrated to verify the effectiveness of the proposed method.


Introduction
As one of the most important problems in sensor array signal processing, direction-of-arrival (DOA) estimation plays an important role in many applications such as radar, sonar, communication system, wireless sensor networks and some other fields. Therefore, it has drawn many researchers' attention and plenty of methods [1] have been proposed to deal with this issue in the past few years, such as MUSIC [2], ESPRIT [3], etc. Most of these methods made the assumption that the sources are located relatively far from the array so that the waves emitted by the sources can be considered plane waves. However, when the sources are located close to the array, the wavefront should be characterized by both the DOA and range. Thus, existing methods based on the far-field assumption are not applicable to this situation.
In the past decades, numerous DOA estimation methods were presented in the near field, such as the two-dimensional MUSIC method [4], maximum likelihood method [4,5], the weighted linear prediction method [6] and higher-order-based methods [7][8][9][10][11]. However, these methods either require multidimensional search or suffer poor resolution from heavy aperture loss. Recently, some DOA estimation methods [12][13][14] based on sparse signal recovery have been developed for near-field sources, exhibiting several advantages compared to other methods. Nevertheless, these sparse-signal-recovery-based methods greatly rely on properly choosing the regularization parameter to reconstruct the sparse signal. Once the regularization parameter is improperly selected, the performance of the sparse-signal-recovery-based methods is greatly influenced. Unfortunately, there is no good way to choose the regularization parameter appropriately and it is still an open issue.

Model of DOA Estimation
Consider a near-field scenario in which K narrowband sources impinge on a symmetric linear array with N = 2u sensors showing Figure 1 (The number of sensor is assumed to be even here, however, our method also works in odd case). The sensors are assumed to be located at the underlying grid of a minimum spacing d. Let x = {x −u , · · · , x −2 , x −1 , x 1 , x 2 , · · · , x u } d be the positions vector of array sensors. The signal received by the ith sensor can be expressed as Equation (1) y i (t) = K ∑ k=1 s k (t)exp(j 2π λ (r ik − r k )) + n i (t), t = 1, 2, · · · , T; i = −u, · · · − 1, 1, · · · , u where Sensors 2016, 16, 1939 3 of 12 denotes the distance between the kth source and the ith sensor, r k represents the range from the kth source to the phase reference, s k (t) is the signal radiated from the kth source, λ denotes the wavelength, n i (t) is the received noise by the ith sensor and T stands for the number of snapshots.
denotes the distance between the th k source and the th i sensor, k r represents the range from the th k source to the phase reference, ( ) k s t is the signal radiated from the th k source,  denotes the wavelength, ( ) i n t is the received noise by the th i sensor and T stands for the number of snapshots. Figure 1. Non-uniform linear array configuration (The number of sensors is assumed to be even).
Using the second-order Taylor expansion of Equation (2), we can approximate 2 ( ) Stacking the measurements of all the sensors in a vector form, we get is the observed signal vector, To make the following derivation simple, some assumptions are made as follows: Using the second-order Taylor expansion of Equation (2), we can approximate τ ik where and Substituting Equation (3) into Equation (1), we have Stacking the measurements of all the sensors in a vector form, we get y(t) = As(t) + n(t), t = 1, 2, · · · , T where y(t) = [y −u (t), y −u+1 (t), · · · , y u (t)] T is the observed signal vector, A = [a(r 1 , θ 1 ), a(r 2 , θ 2 ), · · · , a(r k , θ k )] denotes the array manifold, a(r k , θ k ) = [e jτ −u,k , e jτ −(u−1),k , · · · , e jτ u,k ] T is the so-called steering vector, s(t) = [s −u (t), s −u+1 (t), · · · , s u−1 (t), s u (t)] T and n(t) = [n −u (t), n −u+1 (t), · · · , n u−1 (t), n u (t)] T are the source and noise vectors, respectively. To make the following derivation simple, some assumptions are made as follows: (A1) The sources are non-Gaussian, and mutually uncorrelated. (A2) The noise is additive Gaussian one, either white or coloured, and independent of the sources.

Vitual Array Model by Exploiting Fourth Order Cumulants
The fourth order cumulant of the measured signals is defined as [18] cum(y m (t), y * n (t), y p (t), y * q (t)) = E y m (t)y * n (t)y p (t)y * Under the above assumptions (A1) and (A2) and by using the properties of cumulants, we have where c 4,sk = cum(s k (t), s * k (t), s k (t), s * k (t)) represents the kurtosis of the kth source signal. Note that the noise term is vanished because the fourth order cumulant of Gaussian noise is zero. Since we are only interested in the DOAs of sources, it is better to remove the ϕ k term and retain the ω k term in Equation (9). Let x n = −x m and x p = −x q , and substituting them to Equation (9), we obtain The right side of Equation (10) behaves like the correlation between the mth and qth sensor output of a virtual far-field array, which is just the difference co-array of the physical array.

The DOA Estimation Method
If a uniform linear array (ULA) is exploited, we can easily construct a matrix C 1 , the (m, q)th entry of which can be given by According to [10], C 1 can be represented in a compact form as T and Λ = diag {c 4,s1 , c 4,s2 , · · · , c 4,sK }. Note that the matrix C 1 behaves like the covariance matrix of the received signals by a far-field array, where the array manifold and steering vector are given by the matrix B and the vector b(ω k ). Then, the conventional MUSIC method [2] can be utilized to find the DOAs of the sources. However, with an array of N elements the maximum number of sources that the method can detect is N − 1, hence it is not applicable to the case of underdetermined DOA estimation. Note that according to Equation (10), the difference co-array of a physical array is acquired. Therefore, non-ULA can be employed so as to further increase the resolution and the number of sources to be detected.

Nested Array
In literature [17], a new array geometry called nested array was proposed to locate far-field sources. Through the use of nested array, up to O(N 2 ) far-field sources with N sensors can be detected while obtaining relatively high resolution. Nested array [17] is basically a concatenation of two ULAs: inner and outer where the inner ULA has N 1 elements with spacing d and the outer ULA has N 2 elements with spacing d such that d = (N 1 + 1)d. More precisely, it is a linear array with sensor location given by S inner = {md, m = 1, 2, · · · , N 1 } and S outer = {n(N 1 + 1)d, n = 1, 2, · · · , N 2 }. It has been shown that the difference co-array of nested array is a filled ULA with 2N 2 (N 1 + 1) − 1 elements [17], which implies that the degrees of freedom (DOF) is 2N 2 (N 1 + 1) − 1. Unfortunately, we cannot exploit nested array directly in the near field since the derivation of Equation (10) is based on a symmetric array.

Symmetric Nested Array
In order to take advantage of nested array in the near-field, a simple approach is to locate two identical nested arrays symmetrically. For convenience, we call it a symmetric nested array (SNA). The sensor location of a SNA is given by From the property of nested array, the difference co-array of symmetric nested array with 2(N 1 + N 2 ) elements contains a ULA with 2N 2 (N 1 + 1) − 1 elements, where N 1 and N 2 represent the number of sensors of the inner and outer ULA in each nested array, respectively.

Compressed Symmetric Nested Array
This subsection introduces a new array geometry which is based on nested array. Note that a nested array consists of two ULAs, where the inner ULA and the first elements of the outer ULA also constitute a ULA. Since ULA is always symmetric, we construct the new symmetric nested array as follows.
First, it is assumed that the inner and outer ULA in a nested array have N 1 and N 2 elements, respectively. Since the inner ULA and the first element of the outer ULA also compose of a new ULA with N 1 + 1 elements, we take the centroid of the new ULA as the phase reference. Then a ULA with N 2 − 1 elements is added to the left of the nested array to make the non-ULA symmetrically. Hence, the final non-ULA has N = N 1 + 2N 2 − 1 elements. Compared with a symmetric nested array, fewer sensors are used in the final non-ULA while achieving the same DOF, so we call it a compressed symmetric nested array (CSNA). More specifically, the sensor location of a CSNA is given by the union of the sets To make the definition clearly, an example of a CSNA with 8 elements is depicted in Figure 2. Apparently a CSNA contains a nested array, so it can achieve the same DOF as the contained nested array. The difference co-array of the two arrays contains the same ULA.  It has been shown that a nested array can attain 2N 2 (N 1 + 1) − 1 DOF in the difference co-array using only N 1 + N 2 elements [17]. Similarly, we can find a systematic way to increase the DOF of the difference co-array in the near-field. The distribution of sensors can be further optimized by finding N 1 and N 2 that maximize the total DOF 2N 2 (N 1 + 1) − 1, under the constraint of fixed total number of sensors, i.e., N = N 1 + 2N 2 − 1. The solution can be verified as seen in Table 1.

DOA Estimation Using CSNA in the Near-Field
As Equation (10) shows, the correlation of the difference co-array of the physical array is acquired. More specifically, if the physical array sensor locations are spatially distributed in such a way that the elements in the set of differences x m − x q , −u ≤ m, q ≤ u represent every integer from 0 to N v , where N v ≤ N(N − 1)/2 is an integer, then we have N v + 1 autocorrelation lags given by Interestingly, these autocorrelation lags are identical to those corresponding to a uniform linear array of N v + 1 elements for the same source scene in the far-field. When a CSNA with N = N 1 + 2N 2 − 1 elements is employed, it can be verified that its difference co-array contains the same ULA with 2N 2 (N 1 + 1) − 1 elements as that of the corresponding nested array, which implies that N v = DOF−1 2 = N 2 (N 1 + 1) − 1 sources can be identified using a CSNA with N = N 1 + 2N 2 − 1 elements.
Then, we can construct a matrix: The matrix C has the same form as the conventional covariance matrix of the output of the ULA with N v + 1 elements whose steering vector can be represented as a(θ k ) = [1, e j2ω k , · · · , e j2N v ω k ]. Then, the conventional MUSIC method can be used to estimate the DOAs. Let U n denote the noise space spanned by the eigenvectors corresponding to the small eigenvalues. The spatial spectrum can be defined as As a result, the DOAs can be obtained by finding the first K peaks of P(θ). Given the measured data {y(t), t = 1, 2, · · · , T}, the proposed DOA estimation method using CSNA can be summarized as the following steps.
When SNA is exploited in Equation (12), we can also construct a virtual covariance matrix like Equation (13) and use the conventional MUSIC method to estimate the DOAs. A subspace-based DOA estimation method using a SNA with N = 2(N 1 + N 2 ) elements can detect N 2 (N 1 + 1) − 1 sources. It is obvious that the method based on SNA can identify more sources than sensors when N 1 ≥ 3 and N 2 ≥ 3.

2.
Since only one half of the difference co-array is employed in our method, the proposed method can detect N 2 (N 1 + 1) − 1 sources with N = N 1 + 2N 2 − 1 elements. Compared with the SNA, fewer sensors are acquired for CSNA to detect the same number of sources.

3.
Regarding the computational complexity of the proposed method, the main cost is in calculating cumulants and eigenvalue decomposition (EVD) of matrix C. Calculation of cumulants and EVD of matrix C requires O(TN 4 ) and O((N v + 1) 3 ), respectively. When different array geometry is utilized, the dimension of matrix C is also different. Without loss of generality, we assume the number of sensors N = 4k, where k is an integer. The value of N v + 1 is presented in Table 2 for ULA, SNA and CSNA respectively. . It is also higher than that of the method in [10], where the main cost is O(TN 4 ) + O(N 3 ). However, the advantages of the proposed method based on the CSNA include high resolution and the ability to detect more sources than sensors.

Simulation Results
In this section, we provide some numerical experiments to illustrate the superior performance of the proposed method based on CSNA. In the following figures, the ULA, S-nested array and CS-nested array stand for the method in [10], and the proposed method using the SNA and CSNA, respectively. First of all, an example will be given to demonstrate the underdetermined ability of the proposed method using the CSNA. Then we make a comparison of the three methods in terms of computational cost, RMSE and resolution ability. The source signals and noise are modeled as an exponential process and Gaussian white noise, respectively. The intersensor spacing of the ULA is d = λ 4 , which is equivalent to the underlying grid of minimum spacing of the CSNA and SNA.  Figure 3 plots the spatial spectrum of the proposed method based on the CSNA shown in Figure 2 with SNR = 20dB and T = 4000. It can be seen that the proposed method can resolve the 9 signals clearly with 8 sensors.
respectively. First of all, an example will be given to demonstrate the underdetermined ability of the proposed method using the CSNA. Then we make a comparison of the three methods in terms of computational cost, RMSE and resolution ability. The source signals and noise are modeled as an exponential process and Gaussian white noise, respectively. The intersensor spacing of the ULA is 4 d   , which is equivalent to the underlying grid of minimum spacing of the CSNA and SNA.

Computational Cost
Next, we consider a case in which two uncorrelated sources radiated from  Table 3. By averaging 100 trials, this experiment is implemented in Matlab on windows 7 and a computer with Intel core i5-4590 CPU with 1000 snapshots and a search grid 0.1°. It can be seen that the time consumptions of the three methods are very close in that the main cost of these methods is in calculating cumulants. Besides, it can be also clearly seen that the proposed method based on CSNA requires a little more time than the other two methods, which verifies the theoretical analysis in Section 3. However, by using the CSNA with 8 sensors, the DOF in the proposed method can be increased to

Computational Cost
Next, we consider a case in which two uncorrelated sources radiated from {−11 • , 14λ} and{8 • , 16λ} impinge upon three different linear array with 8 elements, which are the CSNA shown in Figure 2, the SNA with sensor number in each level N 1 = N 2 = 2 and a ULA. The computational time of the above three methods for finding the DOAs are illustrated in Table 3. By averaging 100 trials, this experiment is implemented in Matlab on windows 7 and a computer with Intel core i5-4590 CPU with 1000 snapshots and a search grid 0.1 • . It can be seen that the time consumptions of the three methods are very close in that the main cost of these methods is in calculating cumulants. Besides, it can be also clearly seen that the proposed method based on CSNA requires a little more time than the other two methods, which verifies the theoretical analysis in Section 3. However, by using the CSNA with 8 sensors, the DOF in the proposed method can be increased to N 2 4 + N − 1 = 23 and thus as many as 11 sources, even more than the number of physical sensors, can be resolved. Our proposed method using the SNA can only identify 5 sources because the DOF is only + N 2 = 11 with 8 physical sensors. Furthermore, the proposed method based CSNA can achieve higher resolution than the other two methods, which will be verified in the Section 4.4.

RMSE versus SNR
Now, let us take a comparison of the three methods by evaluating the RMSE of the DOA estimates as a function of SNR. The RMSE is defined as where N mc denotes the number of Monte Carlo trials andθ q,i and θ q,i denote the estimate DOA and the real DOA of the qth signal in the ith trial. Consider two sources located at (−25.4 • , 21λ) and (14.7 • , 13λ) impinge on the three linear arrays with 8 elements used in previous experiment. Figure 4 shows the RMSE versus SNR of the three methods by averaging 200 Monte Carlo trials and calculating 2000 snapshots. It can be seen that the three methods has nearly the same variance in low SNR and the method in [10] shows better performance than the other two methods in moderate and high SNR. Besides, it is easy to discover that the proposed method based on CSNA has lower estimation error compared with the method based on SNA in moderate and high SNR.   impinge on the three linear arrays with 8 elements used in previous experiment. Figure 4 shows the RMSE versus SNR of the three methods by averaging 200 Monte Carlo trials and calculating 2000 snapshots. It can be seen that the three methods has nearly the same variance in low SNR and the method in [10] shows better performance than the other two methods in moderate and high SNR. Besides, it is easy to discover that the proposed method based on CSNA has lower estimation error compared with the method based on SNA in moderate and high SNR.

Resolution Ability versus SNR and the Number of Snapshots
To investigate the resolution ability of the proposed method using CSNA, two closely spaced signals with DOAs

Resolution Ability versus SNR and the Number of Snapshots
To investigate the resolution ability of the proposed method using CSNA, two closely spaced signals with DOAs (−4 • , 4 • ) were impinged on each array. The two sources are considered to be resolved in a trial if both θ 2 − θ 2 and θ 1 − θ 1 are smaller than |θ 2 − θ 1 | /2, where θ 1 and θ 2 denote the true DOAs andθ 1 andθ 2 denote the estimating DOAs. The detection probability of the three methods versus SNR is depicted in Figure 5, where the number of snapshots is T = 2000 and 600 Monte Carlo trials are carried out. As Figure 5 shows, the method using the CSNA gets a higher resolution ability compared with the other two methods, mainly because the CSNA can achieve higher DOF than the others. Figure 6 illustrates the detection probability of the three methods as a function of the number of snapshots, where the parameter are kept the same as before except SNR = 25 dB and DOAs = (−3 • , 3 • ). From Figure 6, it can be clearly observed that, the proposed method based on CSNA outperforms the two other methods due to higher DOF.
higher resolution ability compared with the other two methods, mainly because the CSNA can achieve higher DOF than the others. Figure 6 illustrates the detection probability of the three methods as a function of the number of snapshots, where the parameter are kept the same as before except SNR = 25 dB and DOAs ( 3 ,3 )     .
From Figure 6, it can be clearly observed that, the proposed method based on CSNA outperforms the two other methods due to higher DOF.    and 600 Monte Carlo trials are carried out. As Figure 5 shows, the method using the CSNA gets a higher resolution ability compared with the other two methods, mainly because the CSNA can achieve higher DOF than the others. Figure 6 illustrates the detection probability of the three methods as a function of the number of snapshots, where the parameter are kept the same as before except SNR = 25 dB and DOAs ( 3 ,3 )     .
From Figure 6, it can be clearly observed that, the proposed method based on CSNA outperforms the two other methods due to higher DOF.

Conclusions
In this paper, a novel underdetermined DOA estimation method is proposed based on fourth-order cumulants to achieve superior resolution in the near field. A new array geometry, called a compressed symmetric nested array, is employed so that as many as N 2 (N 1 + 1) − 1 sources can be detected using only N = N 1 + 2N 2 − 1 sensors in the near field. Although our proposed method leads to somewhat higher variance than some other methods due to exploiting parts of the cumulants, it obtains higher resolution. In future work, we will try to utilize all fourth-order cumulants to improve the estimation accuracy.