An improved two-dimensional direction of arrival estimation method of array signal in multipath environment

Abstract

The application provides an improved two-dimensional direction of arrival estimation method of an array signal in a multipath environment, which comprises the following steps: step 1, performing specific sorting on output data of a signal receiving array to obtain non-defective data; step 2, performing low-rank decomposition on the non-defective data in combination with a parallel factor model; and step 3, performing calculation in combination with a least square method criterion and spatial rotation invariance to obtain two-dimensional direction of arrival estimation, and completing improved two-dimensional direction of arrival estimation in a multipath environment. The method provided by the application can provide high-resolution 2D-DOA estimation of a coherent source; the method is suitable for a single fast shot scene; the array element spacing in the method is not limited to half a wavelength; and the method is suitable for a URA array structure and is closer to an actual application scene.

Description

TECHNICAL FIELD

[0001] The present application relates to a method for two-dimensional direction of arrival estimation of array signals, in particular to an improved method for two-dimensional direction of arrival estimation of array signals in multipath environment. BACKGROUND

[0002] Two-dimensional direction-of-arrival (2D-DOA) estimation has been a hot topic in array signal processing for decades. Many algorithms have been proposed for 2D-DOA estimation, such as the algorithms for traditional scalar sensors, like the algorithm of Liu and Mendel (T. H. Liu, J. M. Mendel, “Azimuth and elevation direction finding using arbitrary array geometries,” IEEE Transactions on Signal Processing, 1998, 46(7): 2061-2065.), the algorithm of Wu and Liao (Y. Wu, G. Liao, “A fast algorithm for 2-D direction-of-arrival estimation,” Signal processing, 2003, 83(8): 1827-1831.), and the algorithm of Xi and Liping (N. Xi, L. Liping, “A computationally efficient subspace algorithm for 2-D DOA estimation with L-shaped array,” IEEE signal processing letters, 2014, 21(8): 971-974.). There are also algorithms for new electromagnetic vector sensors (EMVS), like the algorithm of Ren and Ma (S. Ren, X. Ma, “2-D unitary ESPRIT-like direction-of-arrival (DOA) estimation for coherent signals with a uniform rectangular array,” Sensors, 2013, 13(4): 4272-4288.), and the algorithm of Ahmed and Xiaofei (T. Ahmed, Z. Xiaofei, “Rectangular array of electromagnetic vector sensors: tensor modelling / decomposition and DOA-polarisation estimation,” IET Signal Processing, 2019, 13(7): 689-699.).But in the modern complex communication environment, multipath effect is inevitable, that is, the receiving array will sample a large number of coherent signals. The above algorithms are all for the case of non-coherent sources, so they are often no longer applicable. Considering the actual situation, some algorithms for coherent source direction estimation have been proposed, and the typical solution strategies mainly include spatial smoothing algorithm, such as Chen's (Y.M. Chen, "On spatial smoothing for two-dimensional direction-of-arrival estimation of coherent signals," IEEE Transactions on Signal Processing, 1997, 45(7): 1689-1696.) and Xu and Liu's (Xu Yougen, Liu Zhiwen, "Simultaneous estimation of wave direction and polarization parameters of coherent signal source for electromagnetic vector sensor array: spatial smoothing method," Journal of Communications. 05 (2004): 28-38.); Covariance matrix reconstruction algorithm, such as Chen and Kwong's (F.J. Chen, S. Kwong, "ESPRIT-like two-dimensional DOA estimation for coherent signals" IEEE Transactions on Aerospace and Electronic Systems, 2010, 46(3): 1477-1484.); Polarization smoothing algorithm, such as Rahamim and Tabrikian's (D. Rahamim, J. Tabrikian, "Source localization using vector sensor array in a multipath environment," IEEE Transactions on Signal Processing, 2004, 52(11): 3096-3103.). However, the above algorithms also have corresponding defects, including the need to perform eigenvalue decomposition on high-dimensional matrices, which brings heavy computational burden, and the estimation result is affected by noise. SUMMARY

[0003] The technical problem solved by the present application is to provide an improved two-dimensional wave direction estimation method for array signals in a multipath environment.

[0004] In order to solve the above technical problems, the present application discloses an improved two-dimensional wave direction estimation method for array signals in a multipath environment, comprising the following steps:

[0005] Step 1, the output data of the signal receiving array is sorted to obtain non-defective data, the specific method is as follows:

[0006] M and N represent the number of rows and columns of electromagnetic vector sensors in the signal receiving array, wherein M and N are integers; λ represents the signal wavelength, and d represents the array element spacing, the position of the receiving array element is set as follows: the array element is located in the set x-o-y plane, and the arrays along the x-axis and y-axis directions are uniform arrays, that is, the array element spacings are the same; K represents the number of far-field signals, K is a positive integer; let θ k , φ k , γ k and η k be the pitch angle, azimuth angle, auxiliary polarization angle and polarization phase difference of the kth signal source relative to the receiving array, (θ k , φ k ) is the two-dimensional direction of arrival of the signal source, wherein k=1, 2, … K; the output data Y of the signal receiving array is represented as:

[0007]

[0008] wherein, represents the Kronecker product, and represents the Khatri-Rao product, and are the spatial response vectors of the x-axis and y-axis receiving arrays to the kth signal, represents an M×1 matrix belonging to the complex domain, represents an N×1 matrix belonging to the complex domain, is the corresponding polarization response vector, represents a 6×1 matrix belonging to the complex domain, s k is the complex envelope of the kth signal; A x =[a x,1 ,a x,2 ,…,a x,K ] represents the spatial response of the array element on the x-axis to the signal, A y =[a y,1 ,a y,2 ,…,a y,K ] represents the spatial response of the array element on the y-axis to the signal, wherein s L represents the signal complex envelope obtained by the Lth snap, represents a K×L matrix belonging to the complex domain, is a noise sample;

[0009] the cosine along the x-axis and y-axis directions is defined as u k =sinθ k cosφk , v k = sin θ k sin φ k The specific forms of the response vectors are as follows

[0010]

[0011]

[0012]

[0013] Wherein, j represents the imaginary unit;

[0014] Definition A y The Kronecker product of B y = A y ⊙ B, wherein B is the polarization response of the signal, and the matrix data array output Y is obtained by vectorization to obtain a new data The specific form is as follows:

[0015] Z = vec (Y) = (S T ⊙ A x ⊙ B y )1 K,1 +N z = y + N z

[0016] Wherein, vec () represents vectorization of data, y is a noise-free data vector, 1 K,1 is a full 1 vector of K rows and 1 column, is the corresponding noise;

[0017] Define an operator Unvec ( ), which works as follows:

[0018] Unvec 6N (y) = B y (S T ⊙ A x ) T

[0019] That is, the column vector containing 6 × M × N × L elements is taken from top to bottom every 6N elements as a column, and sequentially taken as the first column, the second column to the last M × L column, to obtain a new matrix with a dimension of 6N × ML;

[0020] Define S T ⊙ A x = A s The Kronecker product of S T and A x , and the specific form of the new matrix model X is as follows:

[0021]

[0022] In the formula, is the corresponding noise after rearrangement of data, and at this time A S and A y ⊙B are full rank matrices.

[0023] Step 2, combine the parallel factor model to perform low rank decomposition on non-defective data, specifically including:

[0024] Step 2-1, use the parallel factor model to analyze non-defective data, specifically including:

[0025] The noise-free part in the new matrix model X is defined as:

[0026]

[0027] Among them, D n () represents a diagonal matrix taking the nth row data as the diagonal element;

[0028] The element x p,g,n in the pth row and nth column of is rewritten into the form of the sum of three matrix products:

[0029]

[0030] Among them, 1≤p≤6, [B p,k represents the pth row and kth column element of B; [A y ] n,k represents the nth row and kth column element of A y ; [A S ] g,k represents the gth row and kth column element of A S , and the above formula is the parallel factor model of x p,l,n , and is defined as a three-dimensional matrix.

[0031] Step 2-2, perform low rank decomposition, specifically including:

[0032] According to the three-linear alternating least squares algorithm and the uniqueness theorem of the parallel factor model, from the noisy data X, the estimates y and s of B, A and A are obtained, and the relationship between them and the corresponding true value is as follows:

[0033]

[0034]

[0035]

[0036] Wherein, Pi is column fuzzy matrix; Omega1, Omega2, Omega3 are diagonal matrix of scale fuzzy, E1, E2, E3 are corresponding estimated error matrix.

[0037] The obtained estimated result And Column fuzzy does not affect the estimation accuracy, and scale fuzzy is eliminated by normalization.

[0038] The diagonal matrix of scale fuzzy satisfies Omega1Omega2Omega3=I, I is unit matrix.

[0039] Step 3, combining the least square method criterion with spatial rotation invariance to calculate, obtain two-dimensional wave direction estimation, complete the improved two-dimensional wave direction estimation in multipath environment, the specific method is as follows:

[0040] Step 3-1, using the least square method criterion and spatial rotation invariance to obtain corresponding data, the specific method is as follows:

[0041] For A y , the central wave direction is estimated by using the least square method, and the least square method solution w k is obtained, the specific method is as follows:

[0042]

[0043] Wherein, U is selection matrix, The phase of the kth column of A y , the calculation method is as follows:

[0044]

[0045]

[0046] Take the second element value w k (2) of w k , obtain the estimated value of v k : For A S , two submatrices A x1 , A x2 are obtained by using selection matrix, at this time, the estimated value of u k is obtained by using spatial rotation invariance:

[0047]

[0048] Wherein, denotes the pseudo-inverse of a matrix, diag{} k denotes the k-th diagonal value in a diagonal matrix, angle[] denotes the phase value.

[0049] Step 3-2, two-dimensional wave direction estimation is performed, and the specific method is as follows:

[0050]

[0051] The final (θ k ,φ k ) is the two-dimensional wave direction of the source.

[0052] Advantages:

[0053] 1. The method provided by the application can provide high-resolution 2D-DOA estimation of coherent sources.

[0054] 2. The method provided by the application is suitable for single-shot scenarios.

[0055] 3. The element spacing of the method provided by the application is not limited to half a wavelength.

[0056] 4. The method provided by the application is suitable for URA array structures, and is closer to actual application scenarios. DETAILED DESCRIPTION

[0057] The above and / or other aspects of the application will become more apparent by describing in detail the application with reference to the accompanying drawings, as follows.

[0058] Figure 1 is a schematic diagram of the overall process of the application.

[0059] Figure 2 is a schematic diagram of the receiving array of the application.

[0060] Figure 3 is a 2D-DOA estimation scatter plot of the estimator provided by the application.

[0061] Figure 4 is a comparison diagram of the average RMSE of 2D-DOA estimation with SNR variation.

[0062] Figure 5 is a comparison diagram of the average RMSE of 2D-DOA estimation with sampling number variation.

[0063] Figure 6 is a comparison diagram of the average RMSE of 2D-DOA estimation with correlation coefficient variation. DETAILED DESCRIPTION

[0064] The application provides an improved two-dimensional direction-of-arrival (2D-DOA) estimation method of array signals in a multipath environment, that is, an improved 2D-DOA estimation method of parallel factor (PARAFAC) in a multipath environment, wherein the receiving antenna is characterized in that the array stream pattern is composed of a uniform rectangular array (URA), and each array element is a complete electro-magnetic vector sensor (EMVS), as shown in the figure. Figure 2 The array elements are located on a set x-o-y plane, and the arrays along the x-axis and y-axis directions are uniform arrays, that is, the array element spacings are the same.

[0065] The signal model in the multipath environment involved in the application is as follows:

[0066] M and N (M and N are integers) are used to represent the number of rows and columns of the receiving array, λ is used to represent the signal wavelength, d = λ / 2 is used to represent the array element spacing, and the position of the receiving array element is assumed as shown in the figure. Figure 2 K (K is a positive integer) is used to represent the number of far-field signals, and let θ k , φ k , γ k , η k be the elevation angle, azimuth angle, auxiliary polarization angle and polarization phase difference of the kth (k = 1, 2, …, K) signal source relative to the receiving array, respectively, (θ k , φ k ) is also called the 2D-DOA of the signal source. The receiving array output can be represented as:

[0067]

[0068] wherein, represents the Kronecker product, represents the Khatri-Rao product, and are the spatial response vectors of the x-axis and y-axis receiving arrays to the kth signal, is the corresponding polarization response vector, s k is the complex envelope of the kth signal; A x =[a x,1 ,a x,2 ,…,a x,K ], A y =[a y,1 ,ay,2 ,…,a y,K ], is noise sample. Define u k = sin θ k cos φ k , v k = sin θ k sin φ k , the specific form of the response vector is as follows

[0069]

[0070]

[0071]

[0072] s k = [s k,1 , s k,2 , …, s k,L ] T Expression (2d)

[0073] An improved two-dimensional wave direction estimation method of array signal in multipath environment, as shown in Figure 1 , comprising the following steps:

[0074] Step 1, the output data of the receiving array is sorted to obtain non-defective data; that is, the rank recovery method adopted by the present application is as follows:

[0075] Define B y = A y ⊙B, the matrix data in expression (1) is obtained by vectorization to obtain a new set of data The specific form is as follows

[0076] Z = vec (Y) = (S T ⊙A x ⊙B y )1 K,1 +N z =y+N z Expression (3)

[0077] Where, vec () represents vectorization of data, y is a noise-free data vector, 1 K,1 is a full 1 vector of K rows and 1 column, is the corresponding noise. Further, we define an operator Unvec (), which works as follows:

[0078] Unvec 6N (y) = B y (S T⊙A x ) T Expression (4)

[0079] That is, take every 6N elements of the column vector y containing 6 x M x N x L elements as a column from top to bottom, and take them as the first column, the second column to the last M x L column in turn, to obtain a new matrix with dimensions 6N x ML. And it is easy to get (S T ⊙A x ) T is a full rank matrix, which is proved as follows:

[0080] Let c i , c i (i, j ∈ {1, 2, …, K}) be two different columns in (S T ⊙A x ), that is,

[0081]

[0082]

[0083] Then

[0084]

[0085] Because any two different columns of A x are linearly independent, so Thus Further, (S T ⊙A x ) any two columns are linearly independent. Since (S T ⊙A x ) and its transpose matrix (S T ⊙A x ) T have the same rank, (S T ⊙A x ) T is no longer a rank-deficient matrix, thus solving the problem of rank deficiency of array output data covariance caused by signal coherence, and achieving the purpose of recovering the rank. As long as M ≥ K, rank((S T ⊙A x ) T ) = K.

[0086] Definition The corresponding covariance matrix R of X is:

[0087]

[0088] In the formula, R x is the covariance matrix of , I 6Nis a 6N×6N identity matrix.

[0089] Step 2: Combine Parallel Factor (PARAFAC) model analysis with low-rank decomposition. The specific method is as follows:

[0090] For the noise-free part of X, define:

[0091]

[0092] in, The element x in the pth row and nth column of the p,g,n It can be rewritten as the sum of three matrix products:

[0093]

[0094] Among them, [B] p,k represents the p-th row and k-th column element of B; [A y ] n,k Indicates A y The nth row and kth column element of S ] g,k Indicates A S The g-th row and k-th column element of . Formula (11) is x p,l,n The PARAFAC model, and thus It is defined as a three-dimensional matrix, that is, It can be regarded as a slice of the PARAFAC model in a certain direction. Similarly, the slices in the other two directions can be written as:

[0095]

[0096]

[0097] According to Kruskal (Kruskal JB, "Three-way arrays: rank and uniqueness of trilinear decompositions, with application to arithmetic complexity and statistics," Linear algebra and its applications 18.2(1977):95-138.), if B, A y , A S k-rank (if any k columns of H are independent, then the k-rank of H is k H =k) ​​satisfies:

[0098]

[0099] B, A y , A S is unique, i.e. B, A y , A S can be uniquely determined. In this paper, we estimate them by Trilinear Alternating Least Square (TALS).

[0100] According to equation (8), its LS fitting is:

[0101]

[0102] where ‖‖ F denotes Forbenius norm. The LS least square estimate of is:

[0103]

[0104] where denote the previous estimate of A y and B, respectively. Similarly, according to equations (10) and (11), we have:

[0105]

[0106]

[0107] Note that the estimate at this time is noisy, and the estimate of B, A y , A S from noisy data X is related to the corresponding true value as:

[0108]

[0109]

[0110]

[0111] where Π is the column blur matrix; Ω1, Ω2, Ω3 are diagonal matrices of scale blur, and satisfy Ω1Ω2Ω3=I; E1, E2, E3 are the corresponding estimation error matrices. For the obtained estimate, column blur does not affect the estimation accuracy, and scale blur can be eliminated by normalization.

[0112] Step 3, obtain 2D-DOA estimation by combining the least square (LS) criterion and spatial rotation invariance, the specific method is as follows:

[0113] for After the scale ambiguity is removed by normalization, the phase of the kth column is taken as

[0114]

[0115] Thus, the LS criterion is used to obtain

[0116]

[0117] where

[0118] The LS solution for is

[0119]

[0120] The estimate of v k is

[0121]

[0122] where denotes the second element of

[0123] Similarly, after the scale ambiguity is removed by normalization for we define where J x1 = [I M-1 , 0], J x2 = [0, I M-1 ], and I L is the L x L identity matrix, and other analogies. According to the property we can obtain

[0124] A x1 = J1·A S = (I L S T )⊙(J x1 A x )=S T ⊙(J x1 A x ) Expression (24)

[0125] A x2 = J2·A S = (I L S T )⊙(J x2 A x )=S T ⊙(J x2 A x ) Expression (25)

[0126] ​According to Wen (Wen Fangqing, “Generalized spatial smoothing in bistatic EMVS-MIMO radar,” Signal Processing 193 (2022): 108406.):

[0127]

[0128] where, denotes the pseudo-inverse of a matrix, diag{·} k denotes the kth diagonal value in a diagonal matrix, and angle[·] denotes the phase value.

[0129] Finally, the 2D-DOA estimation is

[0130]

[0131] Embodiment:

[0132] To verify the effectiveness of the framework, the Monte Carlo method is used to evaluate the estimation performance. Here, a MxN URA receiving array is set, and each element is a complete co-sited EMVS. The element spacing is d = λ / 2, λ is the wavelength of the electromagnetic signal (the inverse of the frequency), L is the sampling number, and α is the correlation coefficient. It is assumed that K = 3 far-field signals, whose parameters are θ = (40°, 10°, 60°), φ = (-15°, 45°, 20°), γ = (12°, 39°, 63°), and η = (33°, 47°, -21°). In addition, it is assumed that L samples have been collected. The results of each figure in the simulation depend on 300 independent trials. In the simulation, the signal-to-noise ratio (SNR) is defined as Y and N are both data matrices in expression (1). The performance evaluation uses the root mean square error (RMSE) evaluation. It is worth mentioning here that the algorithm simulation uses an auto regression process for modeling the correlation source intensity, and the correlation coefficient matrix is defined as:

[0133]

[0134]

[0135] where α ∈ [0, 1] is the correlation coefficient. This correlation form is closer to the actual situation and has strong practicality.

[0136] First, the scatter plot results of the 2D-DOA estimation of the proposed coherent estimator are Figure 3Given, the left middle of the figure M = 6, N = 6, SNR = 20dB, L = 1, d = λ / 2, α = 1; the right middle of the figure M = 6, N = 6, SNR = 1dB, L = 200, d = λ / 2, α = 1. It can be clearly seen that all parameters are correctly estimated and automatically paired. The results show that the framework scheme can provide a closed-form solution for 2D-DOA estimation of the signal source, and can provide good DOA estimation in single sampling and high noise conditions.

[0137] Secondly, Figure 4 The average RMSE performance of 2D-DOA estimation of the signal source under different signal-to-noise ratios is given in the middle of the figure, where M = 6, N = 6, L = 500, d = λ / 2, α = 1. In order to highlight the reliability of the scheme, the scheme is compared with Spatial-Smoothing algorithm (Xu Yougen, Liu Zhiwen, "Simultaneous estimation of wave direction and polarization parameters of electromagnetic vector sensor array coherent signal source: spatial smoothing method," Journal of Communications. 05 (2004): 28-38.), ESPRIT-Like algorithm (F. Wen, J. Shi, Z. Zhang. "Closed-form estimation algorithm for EMVS-MIMO radar with arbitrary sensor geometry," Signal Process., 2021, 186, 108117.), Polarization-Smoothing algorithm (D. Rahamim, J. Tabrikian, R. Shavit. "Source localization using vector sensor array in a multipath environment," IEEE Transactions on Signal Processing, 2004, 52(11): 3096-3103.) and the Cramer-Rao bound of URA (marked as CRB). It is worth noting that when the signal-to-noise ratio is low (such as -6dB≤SNR≤2dB), the RMSE of the method proposed in the present application is lower than all the compared algorithms, and accurate direction finding estimation can be achieved (other algorithms fail), which shows that the noise suppression of the present application scheme is very strong and has strong practicability.

[0138] In addition, the performance of the present application under three different conditions is compared, in turn (1) different algorithms under different sampling times (M = 6, N = 6, SNR = 18dB, d = λ / 2, α = 1); (2) different signal correlation strengths (M = 6, N = 6, L = 500, SNR = 2dB, d = λ / 2). The results are shown in the figure, respectively, from Figure 5It can be seen that the RMSE of all algorithms will decrease when the number of samples increases. However, the RMSE of the proposed method is lower than that of all the compared algorithms, which indicates that the proposed scheme can provide more accurate estimation performance. From Figure 6 It can be seen that the proposed algorithm is not sensitive to the signal correlation strength, i.e., it can achieve effective direction finding estimation for both non-coherent and coherent sources.

[0139] In a specific implementation, the present application provides a computer storage medium and a corresponding data processing unit, wherein the computer storage medium can store a computer program, and the computer program can run the invention content of the improved two-dimensional direction of arrival estimation method of array signals in multipath environment and part or all steps in each embodiment when executed by the data processing unit. The storage medium can be a magnetic disk, an optical disk, a read-only memory (ROM), a random access memory (RAM), or the like.

[0140] Those skilled in the art can clearly understand that the technical solutions in the embodiments of the present application can be realized by means of a computer program and its corresponding general hardware platform. Based on such understanding, the technical solutions in the embodiments of the present application can be embodied in the form of a computer program, i.e., a software product, which can be stored in a storage medium and includes a plurality of instructions for causing a device (which can be a personal computer, a server, a single-chip microcomputer, a MUU, or a network device) containing a data processing unit to execute the method described in each embodiment or some parts of the embodiments of the present application.

[0141] The present application provides an improved two-dimensional direction of arrival estimation method of array signals in multipath environment. There are many methods and approaches to implement this technical solution, and the above description is only the preferred embodiment of the present application. It should be noted that for ordinary skilled persons in the technical field, some improvements and refinements can be made without departing from the principles of the present application, and these improvements and refinements should also be considered within the protection scope of the present application. The components not explicitly described in the embodiments can be implemented by using existing technology.

Claims

1. An improved two-dimensional direction of arrival estimation method for array signals in a multipath environment, characterized in that: The following steps are involved: Step 1, performing a specific sort on the output data of the signal receiving array to obtain non-rank-deficient data; Step 2: Combine the parallel factor model to perform low-rank decomposition on the non-rank-deficient data; Step 3: Combine the least squares criterion with spatial rotation invariance to calculate and obtain a two-dimensional direction of arrival estimate, thus completing the improved two-dimensional direction of arrival estimation in a multipath environment. The specific method for performing a specific sorting on the output data of the signal receiving array in step 1 is as follows: M and N are integers representing the number of rows and columns of electromagnetic vector sensors in the signal receiving array. λ represents the signal wavelength and d represents the array element spacing. The positions of the receiving array elements are set as follows: the array elements are located on the set xoy plane, and the arrays along the x-axis and y-axis directions are uniform, that is, the array element spacing is the same. K represents the number of far-field signals, and K is a positive integer. Let θ k ,φ k , γ k and η k are the elevation angle, azimuth angle, auxiliary polarization angle and polarization phase difference of the kth source relative to the receiving array, and their values ​​are in the range of [-90°, 90°], (θ k ,φ k ) is the two-dimensional direction of arrival of the signal source, where k = 1, 2, ... K; the output data Y of the signal receiving array is expressed as: in, represents the Kronecker product, ⊙ represents the Catalina product, and are the spatial response vectors of the receiving arrays on the x-axis and y-axis to the k-th signal, represents an M×1 dimensional matrix belonging to the complex field, represents an N×1 dimensional matrix belonging to the complex field, is the corresponding polarization response vector, represents a 6×1 dimensional matrix belonging to the complex field, s k is the complex envelope of the kth signal; A x =[a x，1 , a x，2 ,…,a x，K ], represents the spatial response of the array element on the x-axis to the signal, A y =[a y，1 , a y，2 ,…,a y，K ], represents the spatial response of the array element on the y-axis to the signal, Among them, s L represents the complex envelope of the signal obtained by the L-th snapshot, represents a K×L dimensional matrix belonging to the complex field, is a noise sample; Define the direction cosines along the x-axis and y-axis as u k = sinθ k cosφ k , v k = sinθ k sinφ k , the specific forms of the above response vectors are as follows Wherein, j represents the imaginary unit; Definition A y With B's Katri Rao Ji for B y =A y ⊙B, where B is the polarization response of the array to the signal, and the array output data Y is vectorized to obtain a new set of data The specific form is as follows: Z=vec(Y)=(S T ⊙A x ⊙B y )1 K，1 +N z =y+N z Among them, vec() means vectorizing the data, y is the noise-free data vector, 1 K，1 is a vector of all 1s with K rows and 1 column, is the corresponding noise; Define an operator Unvec(), which works as follows: Unvec 6N (y)=B y (S T ⊙A x ) T That is, take the column vector y containing 6×M×N×L elements from top to bottom and take every 6N elements as a column as the first column, the second column to the last M×L column, to obtain a new matrix with a dimension of 6N×ML; Define S T ⊙A x =A S For S T With A x The Katrila product of , the specific form of the new matrix model X is as follows: Where, is the corresponding noise after the data is rearranged, and at this time A S With A y ⊙B are all full-rank matrices; The low-rank decomposition of non-rank-deficient data using the parallel factor model described in step 2 specifically includes: Step 2-1, use the parallel factor model to analyze non-rank-deficient data; Step 2-2, perform low-rank decomposition; The analysis of non-rank-deficient data using the parallel factor model described in step 2-1 specifically includes: The noise-free part of the new matrix model X Defined as: in, D n () represents a diagonal matrix with the nth row of data as the diagonal elements; Will The element x in row p and column n p，g，n Rewritten as the sum of three matrix products: Where 1≤p≤6, [B] p，k represents the p-th row and k-th column element of B; [A y ] n，k Indicates A y The nth row and kth column element of S ] g，k Indicates A S The g-th row and k-th column element of p，ml，n The parallel factor model of Defined as a three-dimensional matrix; The low-rank decomposition described in step 2-2 specifically includes: According to the trilinear alternating least squares algorithm and the uniqueness theorem of the parallel factor model, we can get the values ​​of B and A from the noisy data X. y and A S Estimates and And the relationship between it and the corresponding true value is as follows: Where ∏ is the column fuzzy matrix; Ω1, Ω2, Ω3 are the diagonal matrices of scale fuzzy, E1, E2, E3 are the corresponding estimation error matrices; The estimated results obtained in step 2-2 and Column ambiguity does not affect the estimation accuracy, and scale ambiguity is eliminated through normalization; The method for obtaining the two-dimensional direction of arrival estimate described in step 3 is as follows: Step 3-1, using the least squares criterion and spatial rotation invariance to obtain the corresponding data; Step 3-2, perform two-dimensional direction of arrival estimation; The corresponding data is obtained by using the least squares criterion and spatial rotation invariance as described in step 3-1. The specific method is as follows: For A y , use the least squares criterion to estimate the central direction of arrival and obtain the least squares solution w k , the specific method is as follows: Among them, U is the selection matrix, A y The phase of the kth column is calculated as follows: Take w k The second element value w k (2) Get v k Estimated value of For A S , using the selectivity matrix to obtain two sub-matrices A x1 , A x2 , then we can use the spatial rotation invariance to get u k Estimated value of: in, represents the pseudo-inverse of the matrix, diag{} k It means taking the kth diagonal value in the diagonal matrix, and angle[] means taking the phase value.

2. The improved method for estimating the two-dimensional direction of arrival of array signals in a multipath environment according to claim 1, characterized in that: The specific method for estimating the direction of arrival described in step 3-2 is as follows: The final (θ k ,φ k ) is the two-dimensional direction of arrival of the signal source.

3. The improved method for estimating the two-dimensional direction of arrival of array signals in a multipath environment according to claim 2, characterized in that: The diagonal matrix of the scale blur described in step 2-2 satisfies Ω1Ω2Ω3=I, where I is the identity matrix.

Images