CAO Chenghu,ZHAO Yongbo,2,*,PANG Xiaojiao,XU Baoqing,and CHEN Sheng

1.National Laboratory of Radar Signal Processing,Xidian University,Xi’an 710071,China;

2.Information Sensing and Understanding,Xidian University,Xi’an 710071,China

Abstract: This paper takes further insight into the sparse geometry which offers a larger array aperture than uniform linear array(ULA)with the same number of physical sensors.An efficient method based on closed-form robust Chinese remainder theorem(CFRCRT) is presented to estimate the direction of arrival (DOA)from their wrapped phase with permissible errors. The proposed algorithm has significantly less computational complexity than the searching method while maintaining similar estimation precision.Furthermore, we combine all phase discrete Fourier transfer(APDFT) and the CFRCRT algorithm to achieve a considerably high DOA estimation precision. Both the theoretical analysis and simulation results demonstrate that the proposed algorithm has a higher estimation precision as well as lower computation complexity.

Keywords: direction of arrival (DOA), wrapped phase, closedform robust Chinese remainder theorem (CFRCRT), all phase discrete Fourier transfer(APDFT).

1.Introduction

Direction of arrival(DOA)estimation of multiple narrowband sources is a major research issue in array signal processing[1–4],wireless communication[5]and so on.The most well-known nonparametric methods including beamforming technique, Capon’s method, and subspace-based method play a fundamental role in the above applications.Beamforming techniques based on the numerous models such as uniform linear arrays (ULA), uniform rectangular arrays(URA),uniform circular arrays(UCA)[6,7]and nested arrays [8], are proposed the earliest to estimate DOA of signal resource. However, there is a fatal drawback in the above methods whose spectrum is restricted by the Rayleigh resolution limit,i.e.,more signal resources in the far field cannot be distinguished successfully when the interval degree of their DOA is less than the beam width of the antenna element.Capon’s method has been a dominant technique owing to resolving sources within a Rayleigh cell,such as multiple signal classification(MUSIC)[9,10]and estimation of signal parameter via rotational invariance technique (ESPRIT) [11]. However, there are three challenges to deal with to achieve supper-resolution. Firstly,the steering vector must be perfectly constructed and the number of the estimated signal resources must be known in advance. That is to say, supper-resolution is subject to the perfect calibrated array. Secondly, some efforts must be made to effectively reduce the computational complexity.Thirdly,a large number of snapshots must be required to achieve supper-resolution.Subsequently, a set of parametric methods based on the maximum likelihood (ML)paradigm[10]are proposed,but an accurate initialization is required to achieve global convergence. Later on, the sparse signal processing methods based on the compressed sensing(CS)theory provide an alternative solution[12,13].These methods have been demonstrated to achieve supperresolution and a high precision with fewer snapshots and a low signal noise rate(SNR).However,for searching optimal matching on an over-complete basis,more extensive computation and complexity do not satisfy its application in the practical engineering. Moreover, the robustness of the algorithm based on CS depends on the selection of the over-complete basis.

Generally speaking,most algorithms above must satisfy the Spatial-Nyquist theorem,i.e.,the spacing between adjacent sensors must be less than half of the wavelength of the incident signal.In reality,this condition is very difficult to be fulfilled because it is impractical for hardware implementations.Even if it satisfies the Spatial-Nyquist theorem, the adjacent sensors also easily create mutual coupling to degrade the DOA estimation performance [14].In fact, the larger an array aperture is, the higher resolution and stronger interference rejection capabilities it will possess. For a ULA, it is worth noting that it obtains a larger aperture via adding sensors.Hence,the non-uniform sparse array becomes a perfect choice to enlarge the array aperture without increasing the number of sensors.Recently some scholars have proposed DOA estimation algorithms for a sparse array, which possess a larger array aperture than the ULA with the same quantity of sensors.In[15],the DOA is estimated from the wrapped phases in the sparse array by robust generalized Chinese remainder theorem (RGCRT). However, no closed-form solution is proposed.Some methods are based on the Chinese remainder theorem (CRT) to determine the DOA from a sparse array,as have been proposed in some studies[16,17].Although these algorithms do not demand that the distance between adjacent sensors must meet the Spatial-Nyquist theorem, they suffer from the computational complexity for the wrapped phase obtained by searching the peak of the MUSIC spectrum,which may seriously restrict the application in practical engineering.Recently,some novel algorithms for DOA estimation in the sparse array were proposed in[18–27].The algorithms utilize the character of sparse array to increase the number of the degree of freedom and precisely estimate the DOA. Nevertheless, they still undergo computation complexity resulting from grid searching.The closed-form robust CRT(CFRCRT) offers another approach to address this problem.It is worth mentioning that the method proposed in this paper can directly obtain the wrapped phase from the all phase DFT(APDFT)spectrum,which can reduce computational complexity effectively without sacrificing its estimation precision.

The APDFT possesses the excellent performance of the initial phase invariance and prevents the spectrum leakage effectively[28], so the APDFT is combined with the improved CFRCRT algorithm to improve DOA estimation precision.

The remaining of this paper is organized as follows.Problem statement on DOA estimation in the sparse array is given in Section 2. DOA estimation in the sparse array model based on the CRT is introduced in detail in Section 3.The APDFT measurement principle is presented in Section 4.The CFRCRT algorithm with APDFT for DOA estimation is developed in Section 5. In Section 6, simulation results follow to prove the validity of the proposed algorithm.Finally,the conclusions are drawn in Section 7.

2.Problem statement

Some assumptions are given to eliminate some interference factors before creating the sparse array model.

(i)The signal source is located in the far field and can be considered as a parallel plane wave when the signal echo reaches the array.

(ii)The signal source is a point source and narrowband signal so that the DOA is unique related to the structure of the array.

(iii) The received noises are additive Gaussian white noises and they are spatially and temporally independent to each other.

Let us consider the sparse array shown in Fig. 1. The array is composed ofLsensors whose inter-sensor spacing isd1,d2,...,dL-1, respectively.The signal source in the far field arrives at every sensor byθ. Generally speaking,the signal reaching theith sensor can be expressed as

whereA,f0andφirepresent signal amplitude,frequency and initial phase,respectively.

Fig.1 Sparse ULA model

At first, suppose that inter-sensor spacingd1=d2=···=dL-1=d. Let us consider the beamformer output from the directionθ.

where the notation (·)Tdenotes the transpose and the notation (·)Hdenotes the Hermitian transpose.a(θ) =is the steering vector.X(t) =s(t)a(θ0)+N(t) is the plane wave signal from directionθ0.N(t) = [n1(t)n2(t)··· nN(t)]Tis the noise vector from the receiver.Generally speaking,the square of the mode of complex gain functionWHa(θ)is defined as

Assuming that the interested direction isθ0, i.e., setW=a(θ0),the antenna pattern is expressed as

Fig.2 shows that the pattern will generate a grating lobe with the interested direction range

Fig. 2 Spatial spectrum with different inter-sensor spacing wavelength ratios(DOA=0°)

3.DOA estimation model based on CRT

whereλis the wavelength of the carrier signal. Note that arbitrary inter-sensor spacing is larger than half of the wavelength. Then, the apparent phase is ambiguous. In other words,the apparent phase is the remainder of the true phase wrapped by 2π,which is formulated as

wherenjis the unknown integer named folding integer.is the measurement value of the phase difference that is named as apparent phase difference.

Fig.3 Sparse array model

From(5)and(6),we can obtain

In order to infer the DOA estimation algorithm based on CRT, both sides of(7)are multiplied by the parameter

whered0anddjare defined as

The parameterCrepresents a non-negative real number which is chosen as needed.Γjis related to the sparse array geometry.Hence,(8)can be simplified as

wherek ∈Z+andis the multiplicative inverse ofΓj0moduloΓjfor 1 ≤j≤L-1.

ProofSinceΓj0andΓj(1 ≤j≤L-1)are pairwisely co-primes, according to Bezout’s theorem, the modular multiplicative inverse ofΓj0moduloΓj(1 ≤j≤L-1),i.e.,can be obtained as follows.

Consider thejth equation of the simultaneous congruences:

The following equation can be easily obtained by multiplying both sides of(13)byΓj,j0and module byΓj.

Combining (12) with (13), the solution ofnj0can be obtained as follows:

Namely the general solution ofnj0can be written as follows:

Substitutenj0into(13)and solve the solution ofnj.

i.e.,(nj0,ni)has the specific form as follows:

Theorem 2Assume allΓj(1 ≤j≤L-1)are pairwisely co-primes and

ProofAccording to Theorem 1,the following formula can be obtained:

Both sides are multiplied byΓ1and consideringfor somej,we have

According to (23) the following formula can be obtained:

The conclusion in(25)is contradictory with the assumption.Hence,holds.From(21),we have

Assume that allΓjfor 1 ≤j≤L-1 are pairwisely co-primes and

if

whereτandMθdenote the maximal remainder error level and the greatest common divisor of the moduli, respectively.Then,we havefor 1 ≤j≤L-1. □

According to Theorem 2,ifεj(1 ≤j≤L-1)stands for the phase measurement error,i.e.,

According to (10), the following formula can be obtained:

The method based on the improved CFRCRT for DOA estimation can be written as follows.

Step 1Construct a residue setaccording tofrom the given erroneous Δφjfor 1 ≤j≤L-1 and the module set{Mj}according to

Step 2Calculate the remainders moduloMθf(wàn)rom the given noisy remainderas follows:

Step 3Calculate the average value offor 1 ≤j≤L-1:

Step 4Calculate the index of the reference remainder:

Step 6Calculate the remainder ofmodulo

whereis the modular multiplicative inverse ofΓj0moduloΓj.

Step 7Calculateat first:

Step 8Then calculatefor 1 ≤j≤L-1:

Step 9Finally calculate

In order to reduce the error,calculate the mean value of

Finally, the DOA of the interested target can be determined as

Compared with the CRT method based on searching,the CFRCRT method has a much simpler form,i.e.,its computational complexity is greatly reduced.Although the DOA of the interested target can be correctly and rapidly determined by the CRCRT, it is not accurate enough in the engineering application because the CFRCRT method is still sensitive to initial phase remainder noises.

4.APDFT

4.1 Traditional DFT algorithm spectral analysis

where the amplitude-frequency characteristic function is given by

and the phase-frequency characteristic function is given by

Fig.4 Traditional 8-point DFT spectrum

From Fig. 4, we can conclude that the accurate frequency and initial phase values of the sequence can be obtained only whenβis an integer number.The phenomenon above is caused by the spectrum leakage and picket-fence effect because the truncated sequence is done periodic extension by applying the traditional DFT. The APDFT method proposed in this paper has the excellent characteristics of the initial phase with unchanged and restraining spectrum leakage. The principle of APDFT will be introduced in detail in Section 4.2.

4.2 APDFT algorithm spectrum analysis

Consider a finite length sequencex(n) (-N+1 ≤n≤N -1) and a given window functionw1(n) (0 ≤n≤N -1).xm(n)can be obtained by multiplyingx(n)with the shiftedw1(n+m)that indicates shiftingw1(n)bymto the left:

Then,define a periodic sequenceas

Given another window functionw2(n)(0 ≤n≤N-1)and multiplywithw2(n):

i.e.,

wherem= 0,1,2,...,N -1,yrepresentsNdifferent sequences overlapped with each other.The window functionw1(n)andw2(n)mentioned above such as rectangle,hamming,triangle can be chosen as needed.The all phase signalxap(n)can be obtained by superimposingym(n).

From(47),we can observe that every sample fromx(0)tox(N -1)are all included in theNdifferent sequences.For a discrete sequence, sampled points represent phase information, then the phase ofxap(n) is the same. Just like the analysis above,theN-point DFT of the sequencexap(n) will indicate that the phase is a constant number that is equal to the initial phase.

Transfer the sequencexap(n) from the time domain to the frequency domain using the traditional DFT.

where the amplitude-frequency characteristic function is given by

and the phase-frequency characteristic function is given by

Comparing Fig.4 with Fig.5,we can get the following conclusions.

Fig.5 All phase 8-point DFT spectrum

The APDFT possesses excellent property of preventing spectrum leakage because the decreasing velocity of the sidelobe of the APDFT is quadratic times as quick as the traditional DFT.

The APDFT possesses the property of the initial phase invariability because the phase spectrum is independent with the pointk.

Obviously, the APDFT can improve the initial phase estimation precision. Hence, this paper will combine the APDFT with the CFRCRT algorithm to effectively improve the precision of the DOA estimation of the interested target.

According to the discussion above, the APDFT needs three timesNpoint-FFT,i.e.,the computational complexity of the algorithm based on the APDFT isO(Nlog2N).

5.Scheme of DOA estimation based on improved CFRCRT with APDFT

From Fig.6, the scheme of the DOA estimation based on the improved CFRCRT with the APDFT is given by the following.

Fig.6 Scheme of method based on improved CFRCRT for DOA estimation

Step 1Determined0andCas needed.

Step 2Implement the signal in every channel by applying the APDFT and the extract initial phaseφj(1 ≤j≤L)from the spectrum.

Step 3Substitute the module set{mj,1 ≤j≤L-1}and the remainder setinto the improved CFRCRT algorithm.

Step 4Determine the DOAθaccording to(39).

6.Simulation and performance analysis

In this section, some numerical simulations are presented to validate the DOA estimation performance of the proposed algorithm. To compare the estimation performance of the proposed algorithm,the classic and latest DOA estimation algorithms based on the sparse array including the RGCRT proposed in[16],the exploiting sparse array motion(ESAM)algorithm proposed in[18],and compresses sparse array(CSA)algorithm proposed in[19]are utilized as references. In order to ensure the comparability of all tested algorithms,the ULA is selected as the receive array model with the transmit wavelength being 0.3 m.In the following simulation,assume that there are two uncorrelated targets located at (θ1,θ2) = [30°,60°], respectively. Unless noted otherwise, the additive noise follows the complex additional Gaussian random process with the mean zero and the variance.In this paper,we define SNR as

wheredenotes the signal mean power.We use the root mean square error (RMSE) to assess the performance of the tested algorithms.The RMSE is defined as follows:

wheredenotes the real DOA of thejth signal source andis the estimated DOA in theith Monte Carlo trail for thejth signal source.

6.1 Measurement accuracy performance

We firstly inspect the measurement precision performance versus different snapshot numbers for aforementioned four algorithms. The element number of the array isL= 3.The adjacent element spacing is set asd1= 0.6 m,d2=0.4 m,respectively.It is obvious that the spacing between adjacent sensors is larger than half of the wavelength.The 64-points APDFT is adopted to capture the wrapped initial phases of two signal sources.The SNR is set as 10 dB.The input snapshot number varies between 0 and 1 200.

Fig.7 obviously indicates the RMSE gradually declines as the number of snapshots increases. Furthermore, it is sufficiently shown from Fig.8 that the proposed algorithm is superior over the other three algorithms for measurement precision of the DOA since the APDFT effectively and accurately extracts the initial phases of each signal source from the finite sampled data.In addition,the measurement precision can achieve 10-3level when the number of snapshots is equal to 200.Hence,the measurement precision can attain the demand of practical engineering requirement.

Fig.7 RMSE versus the number of snapshots

Secondly, the experiment is conducted to validate the measurement precision of the proposed algorithm versus different SNRs where the number of snapshots is set by 200.The input SNR varies from –10 dB to 10 dB.

Fig.8 evidently exhibits the RMSE declines as the SNR increases.More importantly,it is clearly seen from Fig.8 that the proposed algorithm is slightly better than the other three algorithms when SNR is between –6 dB and 10 dB because the APDFT in the proposed algorithm possesses the excellent performance of the initial phase invariance and prevents the spectrum leakage effectively.In addition,since the RGCRT algorithm measures the phase by grid searching, the estimation precision of the RGCRT algorithm in this paper is slightly better than that of the proposed algorithm when the SNR is less than –6 dB.However,the CSA algorithm and ESAM algorithm depend on some features of the spares array to estimate the DOA whereas they are always not accurate, so the DOA estimation precision is slightly worse.

Fig.8 RMSE versus SNR

As is discussed above, the aforementioned four algorithms have nearly identical performance for estimation precision in the same simulation condition.

6.2 Computational complexity

In this simulation, we consider the computational complexity of the proposed algorithm.The cumulative time is recorded through 50 Monte Carlo trails.It is assumed that the sparse array collects sampling data by 200 snapshots in each Monte Carlo trail and the SNR is set as 10 dB.

It is clearly indicated from Fig.9 that the computational complexity of the proposed algorithm is evidently less than the other three algorithms due to the closed-form analytical solution. The aforementioned three tested algorithms all adopt the grid searching method to estimate the DOA.Hence,the three algorithms cost more time.

Fig.9 Computational complexity of four algorithms

We can draw an important conclusion that the proposed algorithm has nearly identical performance in estimation precision, whereas the computational complexity of the proposed algorithm is much lower than that of the other three tested algorithms.

6.3 Sparse performance

In the following simulation, all of the simulations adopt the proposed algorithm with the 64-points APDFT and 500 snapshots. We first consider the influence of the element number of the sparse arrays on DOA estimation precision.In this simulation,set the adjacent-sensor spacing as[2.4 m, 2.2 m], [2.6 m, 2.4 m, 2.2 m] corresponding toΓ1=[11,12],Γ2=[11,12,13],respectively.

It can be obviously seen from Fig. 10 that the bigger the element number is, the more accurate the DOA estimation will be.Because more information can be provided by more sensors existing in the array, the algorithm can achieve a much higher precision.In addition,the common sense has been proved that the larger the aperture on account of more numbers of the sensors with the same spacing is,the higher the estimation precision will be.

Fig.10 Relation between element number and estimation precision versus SNR(N =64)

Finally, we consider the relation between estimation precision and the sparsity of the array. Set the adjacentsensor spacing as[1.2 m,1 m],[2 m,1.8 m],[2.8 m,2.6 m]respectively corresponding toΓ1= [5,6],Γ2= [9,10],Γ3=[13,14].We make 50 independent Monte Carlo trails and find the following results.

From Fig.11,it can be easily found that the proposed algorithm is appropriate for the arbitrary sparse array in estimating the DOA.In addition,since the larger the adjacentsensor spacing is, the smaller the algorithm’s tolerance to the error of the phase margin is, the DOA cannot be correctly determined in the sparsest array. It is worth mentioning that although the sparser the array is,the larger the aperture will be,the capability of fault-toleration is a dominant factor for DOA estimation precision in this case.

Fig.11 Relation between the sparsity and estimation precision versus SNR(N =64)

7.Conclusions

In this paper, the improved CFRCRT algorithm with the APDFT is proposed to determine the multi-target DOA through solving the phase ambiguity in the sparse array that provides a larger array aperture than ULA with the same number of sensors.Both the theoretical analysis and simulation results demonstrate that DOA estimation based on the improved CFRCRT algorithm with the APDFT can achieve a much higher precision with maintaining a lower computational complexity in the practical engineering.

The CRT algorithm, which is a powerful tool to solve the ambiguity problem,is proposed in even more articles.However, it is seldom proposed that the estimated parameters, which are reconstructed using the CRT algorithm when the residue sets have some limited errors without any auxiliary condition, are applied to the track-before-detect(TBD)technique[31–35].Hence,it calls for further study that the parameters of multiple targets,which are uniquely determined from the erroneous residue sets by a new robust theorem based on CRT,are applied to TBD to improve the estimated performance of multi-target trajectories.