Direction-of-Arrival Estimation Based on Sparse Recovery with Second-order Statistics

Traditional direction-of-arrival (DOA) estimation techniques perform Nyquist-rate sampling of the received signals and as a result they require high storage. To reduce sampling ratio, we introduce level-crossing (LC) sampling which captures samples whenever the signal crosses predetermined reference levels, and the LC-based analogto-digital converter (LC ADC) has been shown to efficiently sample certain classes of signals. In this paper, we focus on the DOA estimation problem by using second-order statistics based on the LC samplings recording on one sensor, along with the synchronous samplings of the another sensors, a sparse angle space scenario can be found by solving an `1 minimization problem, giving the number of sources and their DOA’s. The experimental results show that our proposed method, when compared with some existing normbased constrained optimization compressive sensing (CS) algorithms, as well as subspace method, improves the DOA estimation performance, while using less samples when compared with Nyquist-rate sampling and reducing sensor activity especially for long time silence signal.


Introduction
Direction-of-arrival (DOA) estimation of propagating plane waves is an extensively studied problem in the field of array signal processing, sensor networks, remote sensing, etc.To determine signal source DOAs using multiple measurements vectors, minimum variance distortionless response (MVDR), and multiple signal classification (MUSIC) algorithms are commonly used [1].By construction, all these traditional DOA estimation methods require Nyquistrate sampling of the received signals, which may result in high storage and bandwidth requirements in many sensing systems.
Compressive Sensing (CS) [2] offers a framework for simultaneous sensing and compression of finite-dimensional vectors enabling a potentially significant reduction in the sampling and computation cost at a limited capability sensing system, which depends on linear dimensionality reduction.And the basis pursuit strategy has been used for formulating the DOA estimation problem as a dictionary selection problem where the dictionary entries are produced by discretizing the angle space and then synthesizing the sensor signals for each discrete angle.Sparseness in angle space implies that only a few of the dictionary entries will be required to match the measurements.
The compressive beamforming (CBF) approach based on 1 minimization [3,4] can only take random projections of the received signals at the sensors and has a model for these as delayed and weighted combinations of multiple signal sources coming for different angles, which is substantially different from MVDR, MUSIC and some other convex optimization approaches based on regularization [5,6] which require Nyquist-rate sampling.Although the CBF method does not require the Nyquist-rate sampling at all sensors, it still needs Nyquist-rate sampling at one sensor, which is called reference sensor (RS).Besides, the CBF method just considers the standard CS discrete framework without referring the hardware realization of compressed sampling for practical application.And some existing universal CS measurement instruments cannot be usually used in practice (at least cannot used for the real-time purpose) because of its time-consuming data collection and the difficulty of physical realization.The random convolution based sampling strategy has been investigated by [7], [8] and [9] to solve this drawback.To avoid hardware realization challenges, we acquire less samples from analog signal directly by using Level-crossing (LC) sampling [10,11,12,13,14,15,16,17], in which an analog input signal is compared with a set of quantization levels (can also be called reference levels) and a sample is produced only when the input analog signal changes enough to cross a level, thus it results in nonuniform sampling and saves dynamic power in the analog-to-digital conversion (ADC) and following DSP for those powered by a very small battery and that involve "bursty" signals with varying activity over time, possibly including long periods of silence.Sampling by LC can mimic the behavior of such input signals.As such, the data collection rate is dictated by the signal itself, rather than Nyquist sampling frequency.One direct benefit of such sam-pling is that it can reduce number of samples without sampling in the non-bursty intervals.Higher instantaneous bandwidth can be offered when LC sampling is performed, and resolution is improved without overall increase in bit rate or power consumption because the significant information in the bursty intervals is sampled.The data transmission rate can be reduced by using LC sampling in communication systems [10,11].
In order to obtain samples efficiently for certain types of input, the levels within the amplitude range of the input need to be appropriately assigned in the ADC.In this work, we assume the input dynamic range is known, and we implement a fixed scheme that uniformly assigns levels in the ADC.Usually, we are interested in reconstructing a single signal from LC samplings.In this paper, we consider the case for DOA estimation using a sensor array, where we try to reconstruct a vector of sources' positions using secondorder statistics based on LC samplings of the RS and synchronous samplings of another sensors.
The paper is organized as follows, In Section 2, we provide the LC sampling scheme and LC ADC architecture.We then introduce the DOA estimation algorithm using secondorder statistics based on multiple measurement vectors and a sparse vector finding in Sec. 3, where we also provide complete algorithmic descriptions and corresponding guaranteed performance analysis results.In Sec. 4, a number of simulation results of our proposed approach is compared with that of the state-of-art CS recovery methods and conventional method.Finally, we conclude with lots of simulation results described on speed signals which are collected using LC sampling of the RS and non-uniform samples of another sensors by keeping synchronous sampling of the RS in Section 5.

Nonuniform Signal Processing Tool 2.1 Level Crossing Sampling Scheme
The Level Crossing sampling scheme (LCSS) is one of the signal-dependent sampling schemes, and it is a better choice for sampling the time-varying signals.For LCSS, a sample on the RS is captured only when the input analog signal x(t) crosses one of the quantization levels which are uniformly spaced by a quantum q, and the samples depend on x(t) variations.The (x n ,t n ) is a sample pair with an amplitude x n and a time t n , x n is exactly equal to one of the quantization levels, and current sampling instant t n can be computed by adding the time elapsed dt n between the t n and t n−1 from the previous instant t n−1 , (1)

LC-based Analog-to-Digital Converter
Usually, the sampling instants are exactly known for conventional sampling, and the ADC number of bits determines the ADC resolution as sample amplitudes are quantization.In this paper, we consider a B-bit (2 B reference levels) LC ADC which equipped with an array of 2 B analog comparators, and the comparators compare with input with corresponding reference levels.Without loss of generality, we assume an amplitude-bound signal x(t) ∈ [-A/2, +A/2] that is T second long, and the LC ADC has 2 B levels with uniform spacing q = A/2 B .Let { 1 , 2 , ..., 2 B } represent all the reference levels used by comparators.The ADC compares the input x(t) with all the reference levels and it will record a level crossing with one of if the following comparison holds for a i : ( That is the LC ADC only records the quantization value Q(s i ) of the true signal s i in the interval [(n − 1)τ, nτ].In order to minimize the consuming power, we can also randomly choose the P of the 2 B comparators are on at any time according to [14], that is, a new set of P reference levels is picked and updated every v seconds (v is a constant), which can be accomplished by a digital circuit that periodically updates the set of on comparators by controlling power supply.Note that, in this paper, we need to add an external circuit to make the LC ADC sampling instants of the RS trigger the sampling instants of another sensors, which ensure the sampling instants on the RS and another sensors are synchronous.

DOA Estimation Method
Assume that we consider a far-field consisting of K sources and a sensor array of (M + 1) sensors with an arbitrary geometry.And we also assume the sensor positions are known and are given by η i = [x i , y i , z i ] T (i = 0, 1, ..., M).Our goal is to determine the DOAs of the signal sources by using the received signals.We assume the RS receives a superposition of the time-domain source signals, x 0 (t) = ∑ K k=1 s k (t) + n 0 (t) and the n 0 (t) is the noise at the RS, then the sensor i observes the time-delayed superimposed source signal, where π k = (θ k , φ k ) is the angle pair consisting of the unknown azimuth and elevation of the source, and ∆ i (π k ) is the relative time delay at the i-th sensor for a source with DOA π k with respect to the RS.Here finding the DOA is equivalent to finding the relative time delay, and the time delay ∆ i in (3) can be determined from geometry: where c is the speed of the propagating wave in the medium.
The source angle pair π k lies in the product of space [0, 2π) θ × [0, π) φ , which must be discretized to form an angle dictionary (sparsity basis).Here, we enumerate a finite set of angles for both azimuth and elevation to generate a set of angle pairs B = {π 1 , π 2 , ..., π N t }, where N t determines the angular resolution.Let b denote the sparsity pattern which selects members of the discretized angle-pair set B. A nonzero positive value at index j of b indicates the presence of a target at the angle pair π j .In particular, assume the angle space is discretized in N t points in all, the sparsity transform matrices {Ψ Ψ Ψ i } M i=1 will be of dimension N s × N t (N s N t ), where N s is the sample number.The Ψ Ψ Ψ i for sensor i can be constructed using proper time shifts of x 0 (t) for each π j in B, and the time shift for sensor i with respect to the RS using (4).
Assume the RS records the signal source using LC sampling as x 0 = [x 0 (t lc1 ), x 0 (t lc2 ), ..., x 0 (t lcT N )], then the sampling data on the sensor i can be described as where t lci (i = 1, ..., T N ) are the LC sample instants, x i (t lc j ) denotes a sample at sample instant t lc j for sensor i, and T N is the number of LC samples.Besides, the j-th column of Ψ Ψ Ψ i corresponding to the time shift of the sampled signal x 0 corresponding to the j-th index of the sparsity pattern vector b, which indicates the proper time shift corresponding to the angle pair π j : The matrix Ψ Ψ Ψ i is the sparsity basis corresponding to all discretized angle pairs B at the i-th sensor.Considering the effect of additive sensor noises, the sparsity pattern vector can be recovered using the Dantzing selector [18] convex optimization problem: where ε is a relaxation variable which makes the true b feasible with high probability, since the formulated problem in ( 7) is a convex optimization problem [19], so we solve it numerically using existing solver [20], and we can obtain a global optimum for the problem (7).
The terms Ψ Ψ Ψ T H and Ψ Ψ Ψ T Ψ Ψ Ψ in the (7) constraint are actually auto-and cross-correlations, respectively.Take two signal sources s 1 (t) and s 2 (t) for example, the recorded LC sampling signal at the RS is Assume the signal amplitudes are equal, and the shifted RS LC sample signal at the i-th sensor is when the assumed DOA is π n , and this time shift of the RS LC sample signal is used to populate the n-th column of the Ψ Ψ Ψ matrix.While the true received sample signal on the i-th sensor is where there are different time shifts for the two signals.For Ψ Ψ Ψ T H we get a column vector whose n-th element is where R 11 is the autocorrelation of signal s 1 (t lc ), R 22 is the autocorrelation of signal s 2 (t lc ), R 12 is the cross-correlation between signal s 1 (t lc ) and s 2 (t lc ).For the matrix Ψ Ψ Ψ T Ψ Ψ Ψ, the element in the n-th row and r-th column is According to the two same assumptions as [3,4]: small cross correlation for signals incoherent assumption, and small autocorrelations for signals decorrelate at small lags assumption.In order to make Ψ Ψ Ψ T H−Ψ Ψ Ψ T Ψ Ψ Ψb small, we should make sure that the large elements in the vector Ψ Ψ Ψ T H are canceled by the large terms in Ψ Ψ Ψ T Ψ Ψ Ψb.According to the assumptions, the two largest elements in Ψ Ψ Ψ T H occur when π n = π 1 and π n = π 2 , because these are two peaks in the autocorrelations R 11 (∆ i (π 1 ), ∆ i (π 1 )) and R 22 (∆ i (π 2 ), ∆ i (π 2 )).When we cancel the element R 11 (∆ i (π 1 ), ∆ i (π 1 )) using the row of Ψ Ψ Ψ T Ψ Ψ Ψ corresponding to π n = π 1 , then the vector b must select the column where π r = π 1 .Likewise, to cancel the element R 22 (∆ i (π 2 ), ∆ i (π 2 )), we use the π n = π 2 row and the π r = π 2 column.And all the other elements will be relatively small according to the above assumptions.
About the constraint parameter ε, it will allow the matching of the two signals at their true DOAs.Then the 1 minimization of the selector vector b will tend to pick the signals whose autocorrelation is large or the larger of the two for different signal amplitudes.Besides, the method can also be extended to the case with K unknown sources at DOAs (θ 1 , φ 1 ), (θ 2 , φ 2 ), ...., (θ K , φ K ) impinging on the array.

Simulation Results
The performance of our proposed approach is evaluated in this section using a linear array of 11 sensors uniformly placed on the x-axis, and the first sensor is selected as RS which is placed to be at the origin.For all the addressed scenarios, a DOA space is discretized into a 1 • angular grid.
To demonstrate our proposed method, we will compare the DOA estimation performance with 1 -SVD [21], CBF algorithm [3] and MUSIC method in respective of DOA estimation performance and sample number.Note that, in the CBF method, assume each sensor takes 20 compressive measurements along with a RS with high-rate sampling frequency (much higher than Nyquist-rate sampling) for constructing the dictionary Ψ Ψ Ψ, and another algorithms use Nyquist sampling of each sensor in our simulation.Besides, to show the advantages of our proposed method, the following evaluation measures are employed to evaluate the signal reconstruction quality of the RS receiving data at LC sampling algorithm: the correlation coefficient (CC), the compressive ratio (CR), and the DOA estimation performance: probability of resolved.CC is used to evaluate the similarity between the original signal and its reconstruction by using the data samples: where x i are re-sampled values of the original signal X = [x 1 , x 2 , ..., x n ], and y i are re-sampled values of the reconstruction signal Y = [y 1 , y 2 , ..., y n ], x and ȳ are sample means of X and Y, s x and s y are the sample standard deviations of X and Y, respectively.The CR between the Nyquist sampling signal and LC sampling signal can be described as where S orig and S lcs represent Nyquist sampling number and LC sampling number respectively.And the probability of resolved is defined as follows, it is said to be resolved, if for any signal with DOA θ k , its estimate θk is such that

Angular Resolution Analysis
First, we will analysis the angular resolution of our proposed method, and the achievable resolution of our approach is evaluated and compared with other tested methods in this subsection, two synthetic speech sources are taken and placed in the far-field of the array, and we set 60 reference levels for the RS ADC, that is we only make 60 comparators work for a 6-bit LC ADC.The two speed signal sources used in our simulation are shown in Fig. 1.The RS signal is the sum of the two source signals, and Fig. 2 gives the part of signal we used for DOA estimation and reconstructed signal using Akima interpolation [22] for the LC samplings.
Considering an SNR of 10 dB, subplots in Fig. 3 show two scenarios with angular separation of ∆θ = [6 • , 18 • ] between two uncorrelated sources, respectively.From Fig. 3, we can see that our approach, MUSIC method outperforms CBF approach and 1 -SVD method in scenarios with closely spaced sources (∆θ = 6 • ), and for a large separation angles (∆θ = 18 • ) situation, all the methods can estimate the sig-     nal DOAs correctly, and thus the simulation curves are overlapped.But, for our sample signal, we just use 896 samples, that is, compared with Nyquist samples (1167 samples using f s = 7350 Hz), the CR is 1.3025 (if we consider the whole signal, the CR is 4.7445).Meanwhile the CC is 0.9893, and the CC using Nyquist sampling is 0.9901, which demonstrate that the LC sampling method can obtain comparable signal reconstruction performance as that of the Nyquist sampling.Besides, our method will show more advantages in respect of using less samples and reducing the sensor activity especially for long time silence signal.Note that, we add the simulation result of our DOA estimation method using Nyquist samplings, which demonstrates our method using LC samplings can reduce sensor activity to save the sensor energy.
Furthermore, we can see in the Fig. 4 that the number of LC samples varies with input, which can be explained when we look at the sample signal in the third subpgraph of Fig. 1.More samples only when the utterance occurs.The LC's adaptive nature prevents it from registering many more samples during quiescent interval where there is no information, and enhances its efficiency.

DOA Estimation Performance
In this subsection, the influence of SNR on the spatial spectrum estimation performance of all tested methods was analyzed in this scenario via probability of resolved criterion in over 500 trials.Fig. 5 shows the probability of resolved of the spatial spectrum estimation at various levels of SNR for three uncorrelated signals θ = [−20 • , 41 • , 47 • ], and the main challenge in this scenario stems from the small separation between two sources θ = [41 • , 47 • ].From Fig. 5, we can see that the MUSIC and our proposed method outperform much better than that of 1 -SVD and CBF approaches, and the performance of our method coincides with that of MUSIC at all SNR levels, while the 1 -SVD was unsuccessful in estimating all DOAs correctly for all values of the SNR, and thus 1 -SVD is not visible in Fig. 5.

Conclusion
In this work, we demonstrate the feasibility of our proposed DOA estimation method by using sparse recovery algorithm with second-order statistics.In our solution, we first obtain the LC samplings on the RS and non-uniform samples of another sensors by keeping synchronous sampling of the RS; then we exploit second-order statistics and the sparsity of the sources in the angle domain, and the obtained sparse pattern by solving an 1 minimization problem determines the number of targets and their corresponding DOAs.We accomplish this by demonstrating that our wireless array scheme is robust against noise in the LC sampling samples, which also can be used to recover the data of the RS.And the fact that all array sensors uses even-based measurements will reduce the amount of data that must be communicated of sensors, which has potential in wireless sensor networks where arrays would be formed from distributed sensors.

Fig. 1 .
Fig. 1.Two input speech signal sources and reference signal.