# Fast DOA Estimation Based on the Transform Domain Weighted Noise Subspace Fitting Algorithm for Generalized Sparse Array

Authors: Chunxi Liu; Zhikun Chen; Dongliang Peng
Journal: International Journal of Antennas and Propagation
 (2021)

---

## Abstract

Compared with uniform arrays, a generalized sparse array (GSA) can obtain larger array aperture because of its larger element spacing, which improves the accuracy of DOA estimation. At present, most DOA estimation algorithms are only suitable for the uniform arrays, while a few DOA estimate algorithms that can be applied to the GSA are unsatisfactory in terms of computational speed and accuracy. To compensate this deficiency, an improved DOA estimation algorithm which can be applied to the GSA is proposed in this paper. First, the received signal model of the GSA is established. Then, a fast DOA estimation method is derived by combining the weighted noise subspace algorithm (WNSF) with the concept of “transform domain” (TD). Theoretical analysis and simulation results show that compared with the traditional multiple signal classification (MUSIC) algorithm and the traditional WNSF algorithm, the proposed algorithm has higher accuracy and lower computational complexity.


---

## Body

## 1. Introduction

With the rapid development of information technology, wireless transmission systems, such as mobile communications, radar, and unmanned aerial vehicles (UAV), have become widely used. As such, antenna technology has become more important [1–5], as it has array antenna direction finding, which is an important branch of antenna technology. Most of the existing arrays are uniform arrays, whose element spacing is less than half wavelength of the incident signal, so as to avoid the ambiguity of angle estimation [6–9]. However, when the number of array elements is limited, the aperture of the uniform array will also be affected, which will lead to lower accuracy of DOA estimation [10, 11]. To solve this problem, scholars have proposed the GSA—a nonequidistant array system—in which the distance between adjacent array elements is more than half of the wavelength of the incident signal [12–15]. Compared with the uniform arrays, the GSA can obtain larger array aperture. Therefore, the estimation accuracy and resolution of the algorithm are effectively improved. The GSA is derived from a uniform array, and it can be obtained by redeploying the array elements of uniform array according to some optimization algorithms (such as simulated annealing algorithm [16, 17] and genetic algorithm [18, 19]). For instance, Chen et al. [20] has proposed an effective method to construct GSA based on the improved genetic algorithm.The most significant feature of GSA is that the spacing between adjacent array elements is unequal. In this case, most conventional DOA estimation algorithms which have strict requirements on array structure, such as ESPRIT algorithm [21–23], will lose their effect. On the other hand, a few algorithms that are less dependent on the array structure, such as the MUSIC algorithm [24–26], do not perform well in the condition of low SNR. At present, researchers generally use the WNSF algorithm for DOA estimation [27]. The WNSF algorithm has no requirement on array structure and can accurately estimate signal DOA under bad conditions such as low SNR or small number of snapshots [28–30]. However, the WNSF algorithm needs spectral peak search. In accurate search or two-dimensional DOA estimation, the total number of points of the spatial spectrum is very huge, which incurs vast amount of computation and hinders the application of the algorithm in practical engineering [31, 32].To reduce the computation, a reduced-dimension search algorithm was proposed in [33, 34], which could reduce the dimensionality of the spectral function and then performed multiple low-dimensional searches. However, this algorithm not only reduced the amount of computation but also reduced the accuracy of parameter estimation. To solve this problem, a concept of transform domain (TD) was proposed [35, 36]. The authors in [35] transformed the received signal into the coarray domain and then iteratively corrected the phase offset between the coarray data and presumed model caused by angle biases according to a closed-form formula. On the other hand, Zhang et al. [36] used spherical Fourier domain to construct the array signal model. However, these methods are only applicable to uniform liner arrays, which means that the method is not universal. Nevertheless, it is a concept worth learning and developing.Motivated by the above facts, in this paper, we combine the concept of “transform domain” with the WNSF algorithm to obtain a fast DOA estimation method which is suitable for the GSA. The procedure of the proposed method is as follows. Firstly, we establish the received data model of the GSA. Secondly, the model is transformed into the TD. Then, the intersection of the noise subspace and its conjugate space in the TD is used to replace the noise subspace in the spectral search function. In this manner, the calculation range of the search function can be reduced by half, and the “rapidity” of the algorithm can be achieved. The theoretical analysis and computer simulation results show that the proposed algorithm has higher accuracy and lower computational complexity than the MUSIC algorithm and traditional WNSF algorithm.


## 2. Signal Model of GSA

Consider a GSA composed ofM elements, as shown in Figure 1. The array elements are arbitrarily distributed in space, and the position coordinates of the array elements are given by xm,ym,zm,m=1,2,…,M.Figure 1
Generalized sparse array.The black dots in Figure1 represent the array elements. The elevation angle θ is defined as the intersection angle between the signal and the positive direction of Z-axis, while the azimuth angle φ is defined as the angle between the projection of the signal onto the XOY plane and the positive direction of the X-axis. Here, the value ranges of θ and φ are both −π/2,π/2. Assume that K narrow-band waves impinge upon the GSA from the elevation angle θk and the azimuth angle φk, the data received by a snapshot of the array can be expressed as(1)xt=Aθ,φst+nt,where xt≜x1t,x2t,⋯,xMtT is the M×1 dimensional array receiving data; st≜s1t,s2t,⋯,sKtT is the signal sampling data; and nt is the additive white Gaussian noise matrix with the same dimension as xt. The array steering matrix Aθ,φ≜aθ1,φ1,aθ2,φ2,…,aθK,φK and the steering vector aθk,φk can be expressed as(2)aθk,φk≜ejβk,1,ejβk,2,…,ejβk,M,where j=−1 is the plural unit. Define λ as the signal wavelength. Then, βk,m is given by(3)βk,m=2πλxmsinθkcosφk+ymsinθksinφk+zmcosθk.The process of estimating the 2D spatial angle of theK signals based on the signal model constructed by (1) is called spatial spectral estimation.


## 3. DOA Estimation Algorithm Based on TD-WNSF

### 3.1. Principle of the Traditional WNSF Algorithm

In practical application, according to the received data fromP snapshots, the spatial correlation is estimated using a time average, and the covariance matrix of the array output is obtained by(4)R^XX=1PXXH,where X is the received data matrix of P snapshots. After performing an eigenvalue decomposition, we can derive the following result:(5)R^XX=USDSUSH+UNDNUNH,where US is the signal subspace composed of eigenvectors corresponding to K larger eigenvalues and UN is the noise subspace composed of eigenvectors corresponding to M−K smaller eigenvalues. Accordingly, we can deduce that the signal subspace and the space formed by the steering vector of the array are the same, and the steering vector space of the array and the noise subspace are orthogonal to each other [37]. This orthogonal relation can be expressed as(6)aHθ,φUN=OM×1,where the symbol ‘H’ represents the conjugate transpose operation, aθ,φ is the M×1 steering vector, and O is the 1×M−K zero vector.Considering that the length of the actual received data matrix is limited and noise is mixed in the actual received data matrix,aθ,φ and UN are not completely orthogonal. In other words, equation (6) is not completely valid. Therefore, the following formula of noise subspace fitting (NSF) is considered [38]:(7)θ,φ=minUNHAθ,φF2=mintrUNHAθ,φAHθ,φUN.Equation (7) can be further extended to a weighted form. The relationship between the noise subspace and the steering vector space of the array is given by(8)UNHAθ,φW1/2=0,where W is the weight. Then, the fitting formula given by equation (7) can be transformed into(9)θ,φ=minUNHAθ,φW1/2F2=mintrUNHAθ,φWAHθ,φUN=mintrWAHθ,φUNUNHAθ,φ.Theoretically, the DOA estimation in equation (9) can be estimated by bringing a 2D search to bear on the ranges of all parameters; however, this is computationally exhaustive.


### 3.2. TD-WNSF Algorithm

When running the traditional WNSF algorithm, an extreme value test should be conducted for each point in the spatial spectrum. With the improvement in search accuracy, the number of spectral points is increased further, leading to a sharp increase in the running time of the algorithm. If we can find a way to compress the extremum search range, then the speed of DOA estimation can also be improved. From this analysis, the following transformation is considered:(10)u=sinθcosφ,v=sinθsinφ,1−u2−v2=cosθ.Combining equations (3) and (10) leads to(11)βk,m=2πλxmuk+ymvk+zm1−uk2−vk2.Moreover, combining equations (2) and (11) leads to(12)au,v=ej2πx1u+y1v+z11−u2−v2ej2πx2u+y2v+z21−u2−v2⋮ej2πxMu+yMv+zM1−u2−v2=ej2πx1−u+y1−v+z11−−u2−−v2ej2πx2−u+y2−v+z21−−u2−−v2⋮ej2πxM−u+yM−v+zM1−−u2−−v2∗=a∗−u,−v,where ∗ represents the conjugate operation. The steering vector aθ,ϕ and the transformation domain steering vector au,v are clearly equivalent in a physical sense. Then, −u,−v in equation (12) can be depicted as a virtual mirror signal source symmetry to u,v. Assume that the weight is given by W=I, where I represents the identity matrix. Moreover, combining equations (6) and (12) leads to(13)aHu,vUN=O,a∗u,vHUN∗=aH−u,−vUN∗=O.In equation (13), the steering vector corresponding to the real signal source u,v is orthogonal to the noise subspace UN, whereas the steering vector corresponding to the virtual mirror signal source −u,−v is orthogonal to the conjugate noise subspace UN∗. If we replace the noise space in the WNSF algorithm with the intersection space of UN and UN∗, as the intersection space is orthogonal to the real steering vector and the virtual steering vector, then the WNSF algorithm can generate an extremum at the real signal source and virtual signal source simultaneously. This characteristic means that the DOA estimation only needs to search half of the u,v domain. Therefore, the purpose of “rapidity” can be achieved. This algorithm of constructing the WNSF spectral function in the TD is called the TD-WNSF algorithm.Solving the intersection space ofUN and UN∗ is essential to constructing the TD-WNSF algorithm. The steps discussed briefly describe the method to find the intersection space of the two subspaces.First, we define a concept called “adjoint solution.” Suppose thatα1,α2,…,αs and β1,β2,…,βt are two groups of vectors of the linear space V and a1,a2,…,as,b1,b2,…,bt is a solution of the equation x1α1+x2α2+…+xsαs=y1β1+y2β2+…+ytβt. Then, a1,a2,⋯,as is called the adjoint solution of a1,a2,…,as,b1,b2,…,bt.Suppose that the noise subspaceUN is represented by α1,α2,…,αM−K and the conjugate noise subspace UN∗ is represented by β1,β2,…,βM−K. Then, the intersection space Uinter of UN and UN∗ can be determined by the adjoint solution of(14)x1α1+x2α2+…+xsαs=y1β1+y2β2+…+ytβt,which is specifically expressed by Definition 1 in this study.Definition 1.
Uinter=UN∩UN∗=k1α1+k2α2+…+ksαs|k1,k2,…,ks is the adjoint solution of the solution of equation (14).
The definition can be proven as follows. Supposeζ=Uinter. From ζ∈UN, we can obtain(15)ζ=k1α1+k2α2+…+ksαs,and from ζ∈UN∗, we can derive(16)ζ=l1β1+l2β2+…+lsβs.
Thus, we have(17)k1α1+k2α2+…+ksαs=l1β1+l2β2+…+lsβs.
Consequently,k1,k2,…,ks is an adjoint solution of equation (14).
By contrast, ifk1,k2,…,ks is the adjoint solution of a solution to that of equation (14), then k1,k2,…,ks,l1,l2,…,ls is the solution of equation (14), which is also the expression of equation (17). Given that the left side of equation (17) belongs to UN, whereas the right side of equation (17) belongs to UN∗, we have(18)ζ=k1α1+k2α2+…+ksαs∈Uinter.
In summary,Uinter=UN∩UN∗=k1α1+k2α2+⋯+ksαs|k1,k2,⋯,ks is the adjoint solution of the solution of equation (14). At this stage, the proof has been completed.
On the basis of on the above algorithm, the TD-WNSF spectrum can be defined as(19)fTD−WNSFu,v=1trWAHu,vUinterUinterHAu,v.
As mentioned, this TD-WNSF spectral function generates extreme values at the real and mirror positions of the signal at the same time. Thus, only a half-domain search in theu,v domain is needed to estimate the signal direction. Then, an accurate search is performed at the vicinity of the spectral peak position and its mirror image position to obtain the precise arrival angle. On the basis of the previous analysis, the proposed algorithm does not entail any array structure requirement and can be applied to the GSA.


### 3.3. Description of Algorithm Step

The implementation steps of the proposed method are summarized as follows:(i)
Step 1: perform eigenvalue decomposition of the array-received data in order to obtain the noise subspace.(ii)
Step 2: calculate the intersection space based on the algorithm of Definition1 and then construct the TD-WNSF spectrum based on equation (19).(iii)
Step 3: search the positive half-spectrum of equation (18) to obtain the estimated value u^i,v^i of the DOA parameter in the u,v domain, where i=1,2,…,K.(iv)
Step 4: replaceaθ,φ in equation (6) with au^i,v^i and for the extreme value test. Among the two elements, the one that satisfies aHu,vUN=0 is the real TD-DOA.(v)
Step 5: substitute the TD-DOA, which was obtained in step 4, into equation (10) to calculate the spatial DOA θ^i,φ^i.(vi)
Step 6: perform an accurate searching in a small area nearθ^i,φ^i. The element that satisfies aHθ,φUN=0 is the real spatial domain DOA.As shown in these steps, the algorithm first implements a rough search process to obtain the TD-DOAu^i,v^i. Then, the inverse trigonometric function is transformed to calculate the rough estimation of the angle θ^i,φ^i. Finally, an accurate search is conducted in the small neighborhood of θ^i,φ^i to obtain the fine estimate θ⌢i,φ⌢i. As such, no angle measurement blurring will occur.


### 3.4. Algorithm Complexity Analysis

The traditional MUSIC algorithm, the traditional WNSF algorithm, and the TD-WNSF algorithm are compared in this study. The array structure used by the three algorithms is shown in Figure1. Consider that K uncorrelated signals impinge upon an GSA of M elements, the number of snapshots is given by L, and the number of search points is given by q. Assume that all the three algorithms perform a rough search, which means that the step size is 1°.For the traditional MUSIC algorithm, the modulusaHθ,φUN2 needs to be calculated for each spectral point, and the dimension of UN is M×M−K. Therefore, the computation of the spectral search of the traditional MUSIC algorithm is qM−KM+1. The computation of the eigenvalue decomposition in the M×M dimensional autocovariance matrix is MK+22. Therefore, the total computation of the traditional MUSIC algorithm is qM−KM+1+MK+22.For the traditional WNSF algorithm, the weightW is usually taken as a unit matrix I in the calculation. The operation used to find the trace traHθ,φUN2 is needed to calculate each spectral point. The dimension of UN is M×M−K, and the computation involved in finding the trace is qM−KM+M. The computation of the matrix’s eigenvalue decomposition is the same as that of the traditional MUSIC algorithm, which is MK+22. Therefore, the total computation of the traditional MUSIC algorithm is qM−KM+M+MK+22.For the TD-WNSF algorithm proposed in this study, the operation used to find the tracetraHθ,φUinter2 is needed to calculate each spectral value point. The intersection space Uinter has a lower dimension denoted by M×M−2K compared with UN, and the search range of the TD-WNSF algorithm is reduced by half. Thus, the computation of the peak search is qM−2KM+1/2, while the computation of the eigenvalue decomposition for the autocovariance matrix in the TD-MUSIC is 2M−1×K−12. Therefore, the total computation of the TD-MUSIC algorithm is qM−2KM+1/2+2M−1K−12.Table1 shows a comparison of the computations of the three algorithms with the number of array elements.Table 1
Computations of the three algorithms.
Array number5678910MUSIC33205136695298481274416000WNSF5480864012712174082282428960TD-WNSF54812702172325445165958We can easily see that, with the increase of the number of array elements, the computational complexity of our algorithm is far lower than the other two algorithms, which reflects the rapidity of our algorithm.


## 3.1. Principle of the Traditional WNSF Algorithm

In practical application, according to the received data fromP snapshots, the spatial correlation is estimated using a time average, and the covariance matrix of the array output is obtained by(4)R^XX=1PXXH,where X is the received data matrix of P snapshots. After performing an eigenvalue decomposition, we can derive the following result:(5)R^XX=USDSUSH+UNDNUNH,where US is the signal subspace composed of eigenvectors corresponding to K larger eigenvalues and UN is the noise subspace composed of eigenvectors corresponding to M−K smaller eigenvalues. Accordingly, we can deduce that the signal subspace and the space formed by the steering vector of the array are the same, and the steering vector space of the array and the noise subspace are orthogonal to each other [37]. This orthogonal relation can be expressed as(6)aHθ,φUN=OM×1,where the symbol ‘H’ represents the conjugate transpose operation, aθ,φ is the M×1 steering vector, and O is the 1×M−K zero vector.Considering that the length of the actual received data matrix is limited and noise is mixed in the actual received data matrix,aθ,φ and UN are not completely orthogonal. In other words, equation (6) is not completely valid. Therefore, the following formula of noise subspace fitting (NSF) is considered [38]:(7)θ,φ=minUNHAθ,φF2=mintrUNHAθ,φAHθ,φUN.Equation (7) can be further extended to a weighted form. The relationship between the noise subspace and the steering vector space of the array is given by(8)UNHAθ,φW1/2=0,where W is the weight. Then, the fitting formula given by equation (7) can be transformed into(9)θ,φ=minUNHAθ,φW1/2F2=mintrUNHAθ,φWAHθ,φUN=mintrWAHθ,φUNUNHAθ,φ.Theoretically, the DOA estimation in equation (9) can be estimated by bringing a 2D search to bear on the ranges of all parameters; however, this is computationally exhaustive.


## 3.2. TD-WNSF Algorithm

When running the traditional WNSF algorithm, an extreme value test should be conducted for each point in the spatial spectrum. With the improvement in search accuracy, the number of spectral points is increased further, leading to a sharp increase in the running time of the algorithm. If we can find a way to compress the extremum search range, then the speed of DOA estimation can also be improved. From this analysis, the following transformation is considered:(10)u=sinθcosφ,v=sinθsinφ,1−u2−v2=cosθ.Combining equations (3) and (10) leads to(11)βk,m=2πλxmuk+ymvk+zm1−uk2−vk2.Moreover, combining equations (2) and (11) leads to(12)au,v=ej2πx1u+y1v+z11−u2−v2ej2πx2u+y2v+z21−u2−v2⋮ej2πxMu+yMv+zM1−u2−v2=ej2πx1−u+y1−v+z11−−u2−−v2ej2πx2−u+y2−v+z21−−u2−−v2⋮ej2πxM−u+yM−v+zM1−−u2−−v2∗=a∗−u,−v,where ∗ represents the conjugate operation. The steering vector aθ,ϕ and the transformation domain steering vector au,v are clearly equivalent in a physical sense. Then, −u,−v in equation (12) can be depicted as a virtual mirror signal source symmetry to u,v. Assume that the weight is given by W=I, where I represents the identity matrix. Moreover, combining equations (6) and (12) leads to(13)aHu,vUN=O,a∗u,vHUN∗=aH−u,−vUN∗=O.In equation (13), the steering vector corresponding to the real signal source u,v is orthogonal to the noise subspace UN, whereas the steering vector corresponding to the virtual mirror signal source −u,−v is orthogonal to the conjugate noise subspace UN∗. If we replace the noise space in the WNSF algorithm with the intersection space of UN and UN∗, as the intersection space is orthogonal to the real steering vector and the virtual steering vector, then the WNSF algorithm can generate an extremum at the real signal source and virtual signal source simultaneously. This characteristic means that the DOA estimation only needs to search half of the u,v domain. Therefore, the purpose of “rapidity” can be achieved. This algorithm of constructing the WNSF spectral function in the TD is called the TD-WNSF algorithm.Solving the intersection space ofUN and UN∗ is essential to constructing the TD-WNSF algorithm. The steps discussed briefly describe the method to find the intersection space of the two subspaces.First, we define a concept called “adjoint solution.” Suppose thatα1,α2,…,αs and β1,β2,…,βt are two groups of vectors of the linear space V and a1,a2,…,as,b1,b2,…,bt is a solution of the equation x1α1+x2α2+…+xsαs=y1β1+y2β2+…+ytβt. Then, a1,a2,⋯,as is called the adjoint solution of a1,a2,…,as,b1,b2,…,bt.Suppose that the noise subspaceUN is represented by α1,α2,…,αM−K and the conjugate noise subspace UN∗ is represented by β1,β2,…,βM−K. Then, the intersection space Uinter of UN and UN∗ can be determined by the adjoint solution of(14)x1α1+x2α2+…+xsαs=y1β1+y2β2+…+ytβt,which is specifically expressed by Definition 1 in this study.Definition 1.
Uinter=UN∩UN∗=k1α1+k2α2+…+ksαs|k1,k2,…,ks is the adjoint solution of the solution of equation (14).
The definition can be proven as follows. Supposeζ=Uinter. From ζ∈UN, we can obtain(15)ζ=k1α1+k2α2+…+ksαs,and from ζ∈UN∗, we can derive(16)ζ=l1β1+l2β2+…+lsβs.
Thus, we have(17)k1α1+k2α2+…+ksαs=l1β1+l2β2+…+lsβs.
Consequently,k1,k2,…,ks is an adjoint solution of equation (14).
By contrast, ifk1,k2,…,ks is the adjoint solution of a solution to that of equation (14), then k1,k2,…,ks,l1,l2,…,ls is the solution of equation (14), which is also the expression of equation (17). Given that the left side of equation (17) belongs to UN, whereas the right side of equation (17) belongs to UN∗, we have(18)ζ=k1α1+k2α2+…+ksαs∈Uinter.
In summary,Uinter=UN∩UN∗=k1α1+k2α2+⋯+ksαs|k1,k2,⋯,ks is the adjoint solution of the solution of equation (14). At this stage, the proof has been completed.
On the basis of on the above algorithm, the TD-WNSF spectrum can be defined as(19)fTD−WNSFu,v=1trWAHu,vUinterUinterHAu,v.
As mentioned, this TD-WNSF spectral function generates extreme values at the real and mirror positions of the signal at the same time. Thus, only a half-domain search in theu,v domain is needed to estimate the signal direction. Then, an accurate search is performed at the vicinity of the spectral peak position and its mirror image position to obtain the precise arrival angle. On the basis of the previous analysis, the proposed algorithm does not entail any array structure requirement and can be applied to the GSA.


## 3.3. Description of Algorithm Step

The implementation steps of the proposed method are summarized as follows:(i)
Step 1: perform eigenvalue decomposition of the array-received data in order to obtain the noise subspace.(ii)
Step 2: calculate the intersection space based on the algorithm of Definition1 and then construct the TD-WNSF spectrum based on equation (19).(iii)
Step 3: search the positive half-spectrum of equation (18) to obtain the estimated value u^i,v^i of the DOA parameter in the u,v domain, where i=1,2,…,K.(iv)
Step 4: replaceaθ,φ in equation (6) with au^i,v^i and for the extreme value test. Among the two elements, the one that satisfies aHu,vUN=0 is the real TD-DOA.(v)
Step 5: substitute the TD-DOA, which was obtained in step 4, into equation (10) to calculate the spatial DOA θ^i,φ^i.(vi)
Step 6: perform an accurate searching in a small area nearθ^i,φ^i. The element that satisfies aHθ,φUN=0 is the real spatial domain DOA.As shown in these steps, the algorithm first implements a rough search process to obtain the TD-DOAu^i,v^i. Then, the inverse trigonometric function is transformed to calculate the rough estimation of the angle θ^i,φ^i. Finally, an accurate search is conducted in the small neighborhood of θ^i,φ^i to obtain the fine estimate θ⌢i,φ⌢i. As such, no angle measurement blurring will occur.


## 3.4. Algorithm Complexity Analysis

The traditional MUSIC algorithm, the traditional WNSF algorithm, and the TD-WNSF algorithm are compared in this study. The array structure used by the three algorithms is shown in Figure1. Consider that K uncorrelated signals impinge upon an GSA of M elements, the number of snapshots is given by L, and the number of search points is given by q. Assume that all the three algorithms perform a rough search, which means that the step size is 1°.For the traditional MUSIC algorithm, the modulusaHθ,φUN2 needs to be calculated for each spectral point, and the dimension of UN is M×M−K. Therefore, the computation of the spectral search of the traditional MUSIC algorithm is qM−KM+1. The computation of the eigenvalue decomposition in the M×M dimensional autocovariance matrix is MK+22. Therefore, the total computation of the traditional MUSIC algorithm is qM−KM+1+MK+22.For the traditional WNSF algorithm, the weightW is usually taken as a unit matrix I in the calculation. The operation used to find the trace traHθ,φUN2 is needed to calculate each spectral point. The dimension of UN is M×M−K, and the computation involved in finding the trace is qM−KM+M. The computation of the matrix’s eigenvalue decomposition is the same as that of the traditional MUSIC algorithm, which is MK+22. Therefore, the total computation of the traditional MUSIC algorithm is qM−KM+M+MK+22.For the TD-WNSF algorithm proposed in this study, the operation used to find the tracetraHθ,φUinter2 is needed to calculate each spectral value point. The intersection space Uinter has a lower dimension denoted by M×M−2K compared with UN, and the search range of the TD-WNSF algorithm is reduced by half. Thus, the computation of the peak search is qM−2KM+1/2, while the computation of the eigenvalue decomposition for the autocovariance matrix in the TD-MUSIC is 2M−1×K−12. Therefore, the total computation of the TD-MUSIC algorithm is qM−2KM+1/2+2M−1K−12.Table1 shows a comparison of the computations of the three algorithms with the number of array elements.Table 1
Computations of the three algorithms.
Array number5678910MUSIC33205136695298481274416000WNSF5480864012712174082282428960TD-WNSF54812702172325445165958We can easily see that, with the increase of the number of array elements, the computational complexity of our algorithm is far lower than the other two algorithms, which reflects the rapidity of our algorithm.


## 4. Simulation and Analysis

Consider a generalized sparse array with 16 array elements, the position distribution of its array elements is shown in Table2.Table 2
The position coordinates of the array.
CoordinateNumberX (m)Y (m)Z (m)10002−0.044870.121770.121313−0.08974−0.243550.242644−0.13461−0.365330.3639550.07681−0.049240.1090660.03194−0.171010.230387−0.01293−0.292790.351708−0.05780−0.414560.4730290.198490.233090.09681100.15362−0.984600.21812110.10875−0.220240.33945120.63880−0.342020.46076130.27529−0.259300.20587140.23043−0.147690.32719150.18556−0.269480.44851160.14069−0.39125−0.39125


### 4.1. Simulation 1: Signal DOA Parameter Estimation

This simulation was implemented to evaluate the effectiveness of the proposed TD-WNSF algorithm. Considering a signal with the following values,θ=60o and φ=20o, impinging upon the GSA. SNR is taken to be 0 dB. The running result of our algorithm is shown in Figure 2. From Figure 2, we can see that the estimation of TD-DOA is (0.17, 0.865). According to equation (10), the spatial angle (60°, 20°) and the TD angle (0.17, 0.865) are equivalent, which illustrates that our algorithm can accurately measure the DOA of the signal.Figure 2
TD-WNSF spectrum.


### 4.2. Simulation 2: Relationship between Algorithm Performance and SNR

This simulation is implemented to compare the DOA estimation performances of the MUSIC, WNSF, and TD-WNSF algorithms. The Monte Carlo number isL = 500, and the SNR shifts from −10 to 20 dB. The root-mean-squared error (RMSE) is defined as(20)RMSEθ=1K∑k=1K1L∑l=1Lθ^k,l−θk2,RMSEφ=1K∑k=1K1L∑l=1Lφ^k,l−φk2,where θ^k,l and φ^k,l represent the estimated values of the Kth signal in the Lth Monte Carlo simulation, respectively. Figure 3 shows the relationship between the RMSE of the azimuth angle and SNR, while Figure 4 shows the relationship between the RMSE of the elevation angle and SNR.Figure 3
RMSE of the azimuth angle.Figure 4
RMSE of the elevation angle.From Figures3 and 4, we can see that the RMSE curves of three methods decrease as SNR increases. However, the RMSE curve of our TD-WNSF algorithm locates below that of MUSIC and WNSF all along, which indicates that the TD-WNSF is superior to MUSIC and RD-MUSIC in a different SNR. This is because our method performs many fine-grained searches in the neighborhood of the TD estimation.


### 4.3. Simulation 3: Relationship between Algorithm Performance and the Angle of Incident Signal

To observe the performance of the new approach more clearly, we compared the algorithm performance with the angle of incident signal. The simulation consists of two parts. The first part entails fixing the elevation angles to 15°, 30°, 45°, 60°, and 75°. The change in azimuth angle is conducted to explore the relationship between RMSE and azimuth angle. The second part involves fixing the azimuth angles to 15°, 30°, 45°, 60°, and 75°. The change in elevation angle is performed to explore the relationship between RMSE and elevation angle.The simulation result of the first part is shown in Figure5.Figure 5
Relationship between RMSE and azimuth angle at different elevation angles. (a) Fixed elevation angle of 15°. (b) Fixed elevation angle of 30°. (c) Fixed elevation angle of 45°. (d) Fixed elevation angle of 60°. (e) Fixed elevation angle of 75°.
(a)(b)(c)(d)(e)Figure5 shows the relationship between the performance of algorithms and the azimuth angle of the incident signal under five typical elevation angles. The comparison of the five subgraphs indicates that the array has high sensitivity to the signals with azimuth angle between 0° and 40°.The simulation result of the second part is shown in Figure6.Figure 6
Relationship between RMSE and elevation angle at different azimuth angles. (a) Fixed azimuth angle of 15°. (b) Fixed azimuth angle of 30°. (c) Fixed azimuth angle of 45°. (d) Fixed azimuth angle of 60°. (e) Fixed azimuth angle of 75°.
(a)(b)(c)(d)(e)Figure6 shows the relationship between the performance of algorithms and the incident elevation angle of the incident signal under five typical azimuth angles. The comparison of the five subgraphs indicates that the array has high sensitivity to the signals with an elevation angle between 10° and 35°.


## 4.1. Simulation 1: Signal DOA Parameter Estimation

This simulation was implemented to evaluate the effectiveness of the proposed TD-WNSF algorithm. Considering a signal with the following values,θ=60o and φ=20o, impinging upon the GSA. SNR is taken to be 0 dB. The running result of our algorithm is shown in Figure 2. From Figure 2, we can see that the estimation of TD-DOA is (0.17, 0.865). According to equation (10), the spatial angle (60°, 20°) and the TD angle (0.17, 0.865) are equivalent, which illustrates that our algorithm can accurately measure the DOA of the signal.Figure 2
TD-WNSF spectrum.


## 4.2. Simulation 2: Relationship between Algorithm Performance and SNR

This simulation is implemented to compare the DOA estimation performances of the MUSIC, WNSF, and TD-WNSF algorithms. The Monte Carlo number isL = 500, and the SNR shifts from −10 to 20 dB. The root-mean-squared error (RMSE) is defined as(20)RMSEθ=1K∑k=1K1L∑l=1Lθ^k,l−θk2,RMSEφ=1K∑k=1K1L∑l=1Lφ^k,l−φk2,where θ^k,l and φ^k,l represent the estimated values of the Kth signal in the Lth Monte Carlo simulation, respectively. Figure 3 shows the relationship between the RMSE of the azimuth angle and SNR, while Figure 4 shows the relationship between the RMSE of the elevation angle and SNR.Figure 3
RMSE of the azimuth angle.Figure 4
RMSE of the elevation angle.From Figures3 and 4, we can see that the RMSE curves of three methods decrease as SNR increases. However, the RMSE curve of our TD-WNSF algorithm locates below that of MUSIC and WNSF all along, which indicates that the TD-WNSF is superior to MUSIC and RD-MUSIC in a different SNR. This is because our method performs many fine-grained searches in the neighborhood of the TD estimation.


## 4.3. Simulation 3: Relationship between Algorithm Performance and the Angle of Incident Signal

To observe the performance of the new approach more clearly, we compared the algorithm performance with the angle of incident signal. The simulation consists of two parts. The first part entails fixing the elevation angles to 15°, 30°, 45°, 60°, and 75°. The change in azimuth angle is conducted to explore the relationship between RMSE and azimuth angle. The second part involves fixing the azimuth angles to 15°, 30°, 45°, 60°, and 75°. The change in elevation angle is performed to explore the relationship between RMSE and elevation angle.The simulation result of the first part is shown in Figure5.Figure 5
Relationship between RMSE and azimuth angle at different elevation angles. (a) Fixed elevation angle of 15°. (b) Fixed elevation angle of 30°. (c) Fixed elevation angle of 45°. (d) Fixed elevation angle of 60°. (e) Fixed elevation angle of 75°.
(a)(b)(c)(d)(e)Figure5 shows the relationship between the performance of algorithms and the azimuth angle of the incident signal under five typical elevation angles. The comparison of the five subgraphs indicates that the array has high sensitivity to the signals with azimuth angle between 0° and 40°.The simulation result of the second part is shown in Figure6.Figure 6
Relationship between RMSE and elevation angle at different azimuth angles. (a) Fixed azimuth angle of 15°. (b) Fixed azimuth angle of 30°. (c) Fixed azimuth angle of 45°. (d) Fixed azimuth angle of 60°. (e) Fixed azimuth angle of 75°.
(a)(b)(c)(d)(e)Figure6 shows the relationship between the performance of algorithms and the incident elevation angle of the incident signal under five typical azimuth angles. The comparison of the five subgraphs indicates that the array has high sensitivity to the signals with an elevation angle between 10° and 35°.


## 5. Conclusions

In this paper, a novel and efficient TD-WNSF algorithm for estimating the DOA was proposed, which can be applied to the GSA. We first established the received data model of the GSA. And then, we exploited the symmetry of the fitting estimator in the TD to reduce the spectral searching range. Compared with the traditional MUSIC and WNSF algorithms, the computational complexity of the proposed algorithm is significantly reduced. The simulation results indicate that the TD-WNSF algorithm has high accuracy and efficiency.



---

*Source: 9952575-2021-07-23.xml*