Robust incoherent distribution source two-dimensional arrival angle estimation method under L-shaped array

Abstract

The invention discloses a robust incoherent distribution source two-dimensional arrival angle estimation method under an L-shaped array, which comprises the following steps of: establishing a generalized array flow pattern receiving model under the L-shaped array, and the L-shaped array consists of two uniform linear arrays which are respectively positioned on an x axis and a z axis of an xoz plane; sampling data output by the uniform linear arrays on the x axis and the z axis at N moments; constructing aperture expansion array cross covariance vectors corresponding to the uniform linear arrays on the x axis and the z axis by using statistical average operation and a conjugate symmetry strategy; according to the cross covariance vector of the two aperture expansion arrays, all direction angle estimated values taking the x axis as the reference and all direction angle estimated values taking the z axis as the reference are obtained, and all two-dimensional arrival angle estimated values are obtained through pairing; the method has the advantages that the method is not only suitable for uniform and non-uniform noise scenes, but also does not depend on prior information of the number of signal sources, and the practicability and flexibility of the algorithm are improved.

Description

Technical Field

[0001] The invention belongs to the technical field of array signal processing, and in particular relates to a two-dimensional arrival angle estimation method for a robust incoherent distribution source under an L-type array. Background Art

[0002] In cellular wireless communication systems, angle of arrival (AOA) estimation plays an important role. It is the key for base stations to perform downlink beamforming to ensure the performance of cellular wireless communication systems. Although there have been several excellent achievements in the existing angle of arrival estimation, most of them are based on the point source model of line-of-sight single-path propagation. However, the wireless signals in actual cellular wireless communication systems are often transmitted along multiple paths under fast fading channels. In this case, incoherent distributed sources can better characterize the actual scene. However, due to the virtual expansion of the incoherent distributed source array steering matrix, the existing point source model-based angle of arrival estimation algorithm cannot be directly applied, and high-performance angle of arrival estimation faces serious challenges.

[0003] In recent years, with the in-depth development of theoretical research and the increasing demand for arrival angle estimation of incoherently distributed sources, some arrival angle estimation techniques based on incoherently distributed source models have also been proposed. Representative algorithms include the weighted subspace algorithm published by M. Bengtsson et al. of the Royal Institute of Technology in Sweden in the IEEE Trans. Signal Processing journal, the rotational invariant subspace algorithm published by S. Shahbazpanahi et al. of the University of Ontario Institute of Technology in Canada and Y. Liu et al. of Harbin Engineering University in the IEEE Trans. Signal Processing and IEEE Signal Processing Letters journals, the beam compression algorithm published by Z. Zheng et al. of the University of Electronic Science and Technology of China in the IEEE Trans. Veh. Technology journal, the reduced rank subspace algorithm published by H. Chen et al. of Ningbo University and F. Gao et al. of Tsinghua University in the IEEE Sensors Journal, IEEE Trans. Aerosp. Electron. Systems and IEEE Trans. Signal Processing journals, etc.

[0004] The above algorithms are all proposed based on the uniform white noise model, and require prior information on the number of sources to distinguish the signal subspace from the noise subspace. However, in practical applications, due to the non-ideality of hardware and receiving channels, the noise may be unknown non-uniform noise, and the number of sources may be unknown. In this context, the performance of existing algorithms has dropped sharply or even failed. In view of this, providing a two-dimensional arrival angle estimation method for incoherently distributed sources that can overcome the shortcomings of existing algorithms in practical applications has important research significance and practical application value. Summary of the invention

[0005] The present invention aims to solve the limitations of the prior art in the estimation of two-dimensional arrival angle of incoherent distributed sources, especially the problem of dependence on prior information on the number of signal sources in a uniform noise environment. The present invention provides a robust two-dimensional arrival angle estimation method for incoherent distributed sources under an L-type array, which is not only applicable to uniform and non-uniform noise scenarios, but also does not rely on prior information on the number of signal sources, thereby improving the practicality and flexibility of the algorithm.

[0006] The technical solution adopted by the present invention to solve the above technical problem is: a two-dimensional arrival angle estimation method for robust incoherent distribution sources under an L-type array, characterized by the following steps:

[0007] Step 1: Establish a generalized array flow receiving model under an L-type array: assume that K far-field narrowband unrelated and incoherent distributed sources are incident on an L-type array consisting of two uniform linear arrays. The two uniform linear arrays are respectively located on the x-axis and z-axis of the xoz plane. The number of array elements is M, M>K, and the array element spacing is 1 / 2 of the carrier wavelength. The two uniform linear arrays have no shared array elements. The array element closest to the coordinate origin of the xoz plane in the uniform linear array on the x-axis is taken as the first array element on the x-axis, and the array element closest to the coordinate origin of the xoz plane in the uniform linear array on the z-axis is taken as the first array element on the z-axis.

[0008] Step 2: In the generalized array flow receiving model, the data output by the uniform linear array on the x-axis is sampled at N times; similarly, the data output by the uniform linear array on the z-axis is sampled at N times;

[0009] Step 3: Using statistical averaging, perform cross-correlation operations on the N-time sampling data of the uniform linear array output on the x-axis and the N-time sampling data of the first array element output on the z-axis to obtain the cross-covariance vector r xz Similarly, using statistical averaging, the cross-correlation operation is performed on the N-time sampling data of the uniform linear array output on the z-axis and the N-time sampling data of the first array element output on the x-axis to obtain the cross-covariance vector r zx ;

[0010] Step 4: Using the conjugate symmetry strategy, construct r xz The cross-covariance vector of the aperture expansion array is r e ; Similarly, using the conjugate symmetry strategy, construct r zx The cross-covariance vector of the aperture expansion array is r f ;

[0011] Step 5: If the number of signal sources K is known a priori, then according to r e and r f , and using a closed-form solution algorithm based on rotation-invariant subspace technology to obtain K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis;

[0012] If the number of signal sources K is unknown, then according to r e and r f , and use the sparse total least squares algorithm to obtain K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis;

[0013] Step 6: Pair the K direction angle estimates based on the x-axis and the K direction angle estimates based on the z-axis to obtain K pairs of two-dimensional arrival angle estimates.

[0014] In step 2, the sampled data x(t) at the tth moment of the data output by the uniform linear array on the x-axis is approximately expressed as: x(t)≈A(θ)c 0 (t)+A′(θ)c 1 (t)+n x (t), the sampling data z(t) at the tth moment of the data output by the uniform linear array on the z-axis is approximately expressed as: z(t)≈A(β)c 0 (t)+A′(β)c 2 (t)+n z (t); where t = 1, ..., N, θ represents the angular variable of the signal emitted by the incoherent distribution source based on the x-axis, β represents the angular variable of the signal emitted by the incoherent distribution source based on the z-axis, A(θ) represents the steering matrix of the uniform linear array on the x-axis, and A(θ) = [a(θ 1 ), ..., a(θ K )],a(θ k ) is the kth column vector of A(θ), a(θ k ) represents θ k The associated steering vector, θ k represents the direction angle of the signal emitted by the kth incoherent distribution source based on the x-axis, k = 1, ..., K, A(β) represents the steering matrix of the uniform linear array on the z-axis, A(β) = [a(β 1 ), ..., a(β K)],a(β k ) is the kth column vector of A(β), a(β k ) represents β k The associated steering vector, β k represents the direction angle of the signal emitted by the kth incoherent distribution source based on the z-axis, (θ k , β k ) represents the two-dimensional arrival angle of the kth incoherently distributed source, Represents the value of a(θ k ) in the variable θ k The new vector obtained by differentiation is Indicates the value of a(β k ) in the variable β k The new vector obtained by differentiation, c 0 (t) represents the received signal vector on the L-shaped array after the signals emitted by all incoherent distribution sources propagate through all their respective transmission paths at time t, c 0 (t) = [c 1，0 (t),…,c K，0 (t)] T , c k，0 (t) is c 0 The kth element in (t), c k，0 (t) represents the signal emitted by the kth incoherent distribution source at the tth time through L k The received signal on the L-shaped array after propagating through the transmission paths, s k (t) represents the signal emitted by the kth incoherent distribution source at the tth time, γ k，l (t) represents the complex path gain of the signal transmitted by the kth incoherent distribution source at the tth time along the lth transmission path, L k represents the number of transmission paths of the kth incoherent distribution source, the superscript “T” represents the transposition operation, c 1 (t) represents the received signal vector on the uniform linear array on the x-axis after the signals emitted by all incoherent distribution sources at time t propagate through all their respective transmission paths and expand in their respective directional angles based on the x-axis. c 1 (t) = [c 1，1 (t), ..., c K，1 (t)] T , c k，1 (t) is c 1 The kth element in (t), c k，1 (t) represents the signal emitted by the kth incoherent distribution source at the tth time through L k transmission paths and propagates at θ kThe received signal on the uniform linear array on the x-axis after angular expansion in the direction, It represents the deviation between the azimuth angle based on the x-axis when the signal emitted by the kth incoherent distribution source propagates along the lth transmission path and the azimuth angle based on the x-axis when it propagates along the line-of-sight path. The mean is 0 and the variance is less than 3 degrees. It is a random variable that obeys uniform distribution or Gaussian distribution. 2 (t) represents the received signal vector on the uniform linear array on the z-axis after the signals emitted by all incoherent distribution sources at time t propagate through all their respective transmission paths and expand in their respective directional angles based on the z-axis. c 2 (t) = [c 1，2 (t), ..., c K，2 (t)] T , c k，2 (t) is c 2 The kth element in (t), c k，2 (t) represents the signal emitted by the kth incoherent distribution source at the tth time through L k transmission paths and propagates in β k The received signal on the uniform linear array on the z-axis after angular expansion in the direction, It represents the deviation between the azimuth angle based on the z-axis when the signal emitted by the kth incoherent distribution source propagates along the lth transmission path and the azimuth angle based on the z-axis when it propagates along the line-of-sight path. The mean is 0 and the variance is less than 3 degrees. It is a random variable that obeys uniform distribution or Gaussian distribution. x (t) represents the Gaussian noise vector of the uniform linear array output on the x-axis at the tth moment, n z (t) represents the Gaussian noise vector output by the uniform linear array on the z-axis at the tth moment, n x (t) and n z (t) are not related to each other.

[0015] In step 3, Among them, z 1 (t) represents the sampling data of the first array element output on the z-axis at the tth time, x 1 (t) represents the sampled data of the data output by the first array element on the x-axis at the tth moment, the superscript “H” represents the conjugate transpose operation, and p represents the power vector of the received signal on the L-type array after the signals emitted by all incoherent distribution sources at N moments propagate through all their respective transmission paths, and p = [p 1 ,…,p K ] T , pk is the kth element of p, p k The signal emitted by the kth incoherent distribution source at N moments is transmitted by L k The power of the received signal on the L-shaped array after propagating through the transmission path, The superscript “*” indicates a conjugation operation.

[0016] In step 4, Among them, r e and r f The dimensions of are (2M-1)×1, J represents a permutation matrix in which the elements on the anti-diagonal are 1 and the other elements are 0. Indicated by r zx The subvector consisting of the 2nd to Mth elements of Indicated by r zx The subvector consisting of the 2nd to Mth elements of .

[0017] In step 5, according to r e and r f , and the specific process of obtaining K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis using a closed-form solution algorithm based on the rotation-invariant subspace technology is as follows:

[0018] Step 5.1a: Set r e Divide into M sub-vectors with overlapping elements, and transform r e The mth sub-vector obtained by partitioning Expressed as Similarly, r f Divide into M sub-vectors with overlapping elements, and transform r f The mth sub-vector obtained by partitioning Where m = 1, ..., M, Indicated by r e The subvector consisting of the mth element to the M+m-1th element of Indicated by r f The subvector consisting of the mth element to the M+m-1th element of ;

[0019] Step 5.2a: According to r e All sub-vectors obtained by partitioning are constructed r xz The associated M×M matrix Similarly, according to r f All sub-vectors obtained by partitioning are constructed r zx The associated M×M matrix Among them, Λ c0 represents the received signal power matrix with a diagonal structure, Λ c0The dimension is K×K, Λ c0 The element in the kth row and kth column is p k ;

[0020] Step 5.3a: Perform eigenvalue decomposition and get The eigenvectors corresponding to the K largest eigenvalues ​​of , the K eigenvectors constitute the matrix U e , then U e Divide into two sub-matrices and Likewise, Perform eigenvalue decomposition and get The eigenvectors corresponding to the K largest eigenvalues ​​of , the K eigenvectors constitute the matrix U f , then U f Divide into two sub-matrices and Among them, U e and U f The dimensions of are all M×K, Indicated by U e The submatrix consisting of the 1st to M-1th rows of Indicated by U e The submatrix composed of the 2nd to Mth rows of Indicated by U f The submatrix consisting of the 1st to M-1th rows of Indicated by U f The submatrix composed of the 2nd to Mth rows of ;

[0021] Step 5.4a: Calculation Then for H e Perform eigenvalue decomposition to obtain K eigenvalues ​​u 1 , ..., u K ; Then through the formula k'∈[1, K] calculates θ k' , the estimated value of k'∈[1,K] k'∈[1,K]; where superscript represents the pseudo-inverse operation of the matrix, u k' Expressing e The k'th eigenvalue obtained by eigenvalue decomposition, ∠u k' Indicates u k' Phase operation, θ k' represents the k'th direction angle based on the x-axis, Represents θ k' An estimated value of

[0022] Similarly, calculation Then for H fPerform eigenvalue decomposition to obtain K eigenvalues ​​v 1 , …, v K , and then through the formula k'∈[1, K] calculates β k' , the estimated value of k'∈[1,K] Among them, v k' Expressing f The k'th eigenvalue obtained by eigenvalue decomposition, ∠v k' Indicates v k' Phase operation, β k' represents the k'th direction angle based on the z-axis, Represents β k' The estimated value of .

[0023] In step 5, according to r e and r f , and the specific process of using the sparse total least squares algorithm to obtain K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis is:

[0024] Step 5.1b: Set r e Sparse representation is r e =(Ψ(θ)+E θ ) Nθ +e θ , r f Sparse representation is r f =(Ψ(β)+E β ) Nβ +e β ; where Ψ(θ) represents the matrix [JA * (θ), A (2:M) (θ)] T The sparse representation of A (2:M) (θ) represents the submatrix consisting of the 2nd to Mth rows of A(θ), Ψ(β) represents the matrix [JA * (β), A (2:M) (β)] T The sparse representation of A (2:M) (β) represents the submatrix consisting of the 2nd to Mth rows of A(β), E θ represents the θ-related deviation matrix caused by the application of generalized array flow pattern, E β represents the β-related deviation matrix caused by the application of generalized array flow pattern, p Nθ and p Nβ are all K sparse vectors, p Nθ and p Nβ contains both non-zero elements and zero elements, p Nθ The index value of the k'th non-zero element in is θ 1,…,θ K One of the Nβ The index value of the k'th non-zero element in is β 1 , …, β K one of them, k'∈[1,K], e θ Represents r e The deviation vector between the estimated value and the true value, e β Represents r f The deviation vector between the estimated value and the true value of

[0025] Step 5.2b: Construct the sparse total least squares optimization problem P1 associated with θ, described as: P1: Similarly, construct the β-related sparse total least squares optimization problem P2, which is described as: P2: Among them, |||| 2 It means to find the l of the vector 2 norm,|||| F It means to find the Frobenius norm of the matrix, |||| 1 It means to find the l of the vector 1 norm, λ represents the regularization parameter, λ>0, Indicates p Nθ The estimated value of Indicates E θ The estimated value of Indicates p Nβ The estimated value of Indicates E β An estimated value of

[0026] Step 5.3b: Based on the θ-related sparse total least squares optimization problem P1 and the β-related sparse total least squares optimization problem P2, use the alternating descent suboptimal algorithm to complete p through alternating iterations Nθ and p Nβ The specific process is as follows:

[0027] Step 5.3.1b: Let iter represent the number of iterations, and the initial value of iter is 1; let iter max Indicates the maximum number of iterations;

[0028] Step 5.3.2b: Given and Calculate E at the iterth iteration θ Estimated value of and E β Estimated value of Among them, when iter=1 For p Nθ The initial estimate of For pNβ The initial estimate of When iter>1 Indicates that p is the iter-1th iteration Nθ The estimated value of Indicates that p is the iter-1th iteration Nβ An estimated value of

[0029] Step 5.3.3b: Given and Calculate p at the iterth iteration Nθ Estimated value of and p Nβ Estimated value of The calculation process is: construct l 2 -l 1 The norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get The calculation process is: construct l 2 -l 1 The norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get

[0030] Step 5.3.4b: Determine whether iter is equal to iter max If yes, the alternating iteration process ends and the obtained and Corresponding to p Nθ The final estimate of p Nβ The final estimated value; otherwise, let iter = iter + 1, and then return to step 5.3.2b to continue execution; where the "=" in iter = iter + 1 is the assignment symbol;

[0031] Step 5.4b: By finding p Nθ The indices of the non-zero elements in the final estimate of θ determine 1 ,…,θ K Estimated value of p Nθ The index value of the k'th non-zero element in the final estimate of is Represents the k'th direction angle θ based on the x-axis k' Similarly, by finding p Nβ The indices of the nonzero elements in the final estimate of β determine 1 , …, β K Estimated value of p NβThe index value of the k'th non-zero element in the final estimate of is Indicates the k'th direction angle β based on the z axis k' The estimated value of .

[0032] In step 5.3.2b, The solution process is: construct l 2 -l 1 The norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get The solution process is: construct l 2 -l 1 The norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get

[0033] The specific process of step 6 is as follows:

[0034] Step 6.1: Calculate the L-array cross-covariance matrix R xz , Among them, Λ c0 represents the received signal power matrix with a diagonal structure, Λ c0 The dimension is K×K, Λ c0 The element in the kth row and kth column is p k ;

[0035] Step 6.2: R xz Perform column vectorization to obtain vector r v ;

[0036] Step 6.3: According to r v , use the greedy algorithm to unambiguously pair the K direction angle estimates based on the x-axis and the K direction angle estimates based on the z-axis, as follows:

[0037] Step 6.3.1: Construct the complete basis matrix Ξ,

[0038]

[0039] ;in, represents the Kronecker product;

[0040] Step 6.3.2: Calculate |Ξ H r v |, get K K × 1 dimensional output results Among them, || means taking the absolute value operation, Represents the output of the k'th K×1 dimension;

[0041] Step 6.3.3: The maximum value of k' and Complete the pairing process, if Maximum, then and is a pair of two-dimensional arrival angle estimates for the same incoherently distributed source.

[0042] Compared with the prior art, the advantages of the present invention are:

[0043] (1) Two cross-covariance vectors are used to suppress the influence of noise and expansion angle, and a conjugate symmetry strategy is used to expand the virtual array aperture, so that the method of the present invention is not only well applicable to uniform and non-uniform noise, but also has good robustness and high estimation accuracy.

[0044] (2) The closed-form solution algorithm based on the rotationally invariant subspace technology and the sparse total least squares algorithm are used to independently estimate the two-dimensional arrival angle, and then a parameter pairing strategy is adopted. This allows the method of the present invention to achieve unambiguous estimation of the two-dimensional arrival angle under the conditions where the number of signal sources is known or unknown. It has strong practical applicability, few restrictive constraints, and good flexibility. BRIEF DESCRIPTION OF THE DRAWINGS

[0045] Figure 1 Schematic diagram of the implementation process of the method of the present invention;

[0046] Figure 2 Schematic diagram of a generalized array flow receiving model under an L-type array in the method of the present invention;

[0047] Figure 3 It is a simulation result diagram of the method of the present invention and the existing rotationally invariant subspace algorithm and reduced rank subspace algorithm under the condition of known number of signal sources and uniform white noise;

[0048] Figure 4 It is a simulation result diagram of the method of the present invention and the existing rotationally invariant subspace algorithm and reduced rank subspace algorithm under the conditions of unknown number of signal sources and non-uniform noise;

[0049] Figure 5 It is a simulation result diagram of the method of the present invention and the existing rotationally invariant subspace algorithm and reduced rank subspace algorithm under the conditions of unknown number of signal sources, non-uniform noise and different angle expansion variances. DETAILED DESCRIPTION

[0050] The present invention will be further described below with reference to the accompanying drawings and in combination with specific implementations, so that those skilled in the art can implement the invention with reference to the description. The protection scope of the present invention is not limited to the specific implementations.

[0051] The embodiment of the present invention relates to a method for estimating two-dimensional arrival angle of a robust incoherently distributed source under an L-type array, and a schematic diagram of its implementation process is shown in FIG. Figure 1 As shown, the steps are as follows:

[0052] Step 1: Establish a generalized array flow receiving model under L-type array: Figure 2 As shown, it is assumed that K far-field narrowband unrelated incoherent distribution sources are incident on an L-shaped array composed of two uniform linear arrays, K≥1, and K=2 in this embodiment. The two uniform linear arrays are respectively located on the x-axis and z-axis of the xoz plane, and the number of array elements is M, M>K, and M=10 in this embodiment. The array element spacing is 1 / 2 of the carrier wavelength, and the two uniform linear arrays have no shared array elements. The array element closest to the coordinate origin of the xoz plane in the uniform linear array on the x-axis is taken as the first array element on the x-axis, and the array element closest to the coordinate origin of the xoz plane in the uniform linear array on the z-axis is taken as the first array element on the z-axis. Figure 2 In the equation, 1, ..., M represents the number of array elements on the x-axis or z-axis. Figure 2 In the figure, the first array element on the x-axis is the array element numbered 1 on the x-axis, and the distance between it and the coordinate origin of the xoz plane is 1 / 2 of the carrier wavelength. The first array element on the z-axis is the array element numbered 1 on the z-axis, and it is at the coordinate origin of the xoz plane, that is, the distance between it and the coordinate origin of the xoz plane is 0. Of course, in actual processing, the first array element on the x-axis can also be at the coordinate origin of the xoZ plane, and the distance between the first array element on the z-axis and the coordinate origin of the xoz plane is 1 / 2 of the carrier wavelength.

[0053] Step 2: In the generalized array flow receiving model, the data output by the uniform linear array on the x-axis is sampled at N moments; similarly, the data output by the uniform linear array on the z-axis is sampled at N moments.

[0054] Further defined, in step 2, the sampling data x(t) at the tth moment of the data output by the uniform linear array on the x-axis is approximately expressed as: x(t)≈A(θ)c 0 (t)+A′(θ)c 1 (t)+n x (t), the sampling data z(t) at the tth moment of the data output by the uniform linear array on the z-axis is approximately expressed as: z(t)≈A(β)c 0 (t)+A′(β)c 2 (t)+n z(t); wherein N>1, t=1, ..., N, θ represents the angular variable of the signal emitted by the incoherent distribution source based on the x-axis, β represents the angular variable of the signal emitted by the incoherent distribution source based on the z-axis, A(θ) represents the steering matrix of the uniform linear array on the x-axis, and A(θ)=[a(θ 1 ), ..., a(θ K )],a(θ k ) is the kth column vector of A(θ), a(θ k ) represents θ k The associated steering vector, θ k represents the direction angle of the signal emitted by the kth incoherent distribution source based on the x-axis, k=1,...,K,a(θ 1 ) is the first column vector of A(θ), a(θ K ) is the Kth column vector of A(θ), A(β) represents the steering matrix of the uniform linear array on the z-axis, A(β)=[a(β 1 ), ..., a(β K )],a(β k ) is the kth column vector of A(β), a(β k ) represents β k The associated steering vector, β k represents the direction angle of the signal emitted by the kth incoherent distribution source based on the z-axis, a(β 1 ) is the first column vector of A(β), a(β K ) is the Kth column vector of A(β), (θ k , β k ) represents the two-dimensional arrival angle of the kth incoherently distributed source, Represents the value of a(θ k ) in the variable θ k The new vector obtained by differentiation is Indicates the value of a(β k ) in the variable β k The new vector obtained by differentiation, c 0 (t) represents the received signal vector on the L-shaped array after the signals emitted by all incoherent distribution sources propagate through all their respective transmission paths at time t, c 0 (t) = [c 1，0 (t), ..., c K，0 (t)] T , c k，0 (t) is c 0 The kth element in (t), c k，0 (t) represents the signal emitted by the kth incoherent distribution source at the tth time through L kThe received signal on the L-shaped array after propagating through the transmission paths, s k (t) represents the signal emitted by the kth incoherent distribution source at the tth time, γ k，l (t) represents the complex path gain of the signal transmitted by the kth incoherent distribution source at the tth time along the lth transmission path, L k represents the number of transmission paths of the kth incoherent distribution source, the superscript “T” represents the transposition operation, c 1 (t) represents the received signal vector on the uniform linear array on the x-axis after the signals emitted by all incoherent distribution sources at time t propagate through all their respective transmission paths and expand in their respective directional angles based on the x-axis. c 1 (t) = [c 1，1 (t), ..., c K，1 (t)] T , c k，1 (t) is c 1 The kth element in (t), c k，1 (t) represents the signal emitted by the kth incoherent distribution source at the tth time through L k transmission paths and propagates at θ k The received signal on the uniform linear array on the x-axis after angular expansion in the direction, It represents the deviation between the azimuth angle based on the x-axis when the signal emitted by the kth incoherent distribution source propagates along the lth transmission path and the azimuth angle based on the x-axis when it propagates along the line-of-sight path. The mean is 0 and the variance is less than 3 degrees. It is a random variable that obeys uniform distribution or Gaussian distribution. 2 (t) represents the received signal vector on the uniform linear array on the z-axis after the signals emitted by all incoherent distribution sources at time t propagate through all their respective transmission paths and expand in their respective directional angles based on the z-axis. c 2 (t) = [c 1，2 (t), ..., c K，2 (t)] T , c k，2 (t) is c 2 The kth element in (t), c k，2 (t) represents the signal emitted by the kth incoherent distribution source at the tth time through L k transmission paths and propagates in β k The received signal on the uniform linear array on the z-axis after angular expansion in the direction, It represents the deviation between the azimuth angle based on the z-axis when the signal emitted by the kth incoherent distribution source propagates along the lth transmission path and the azimuth angle based on the z-axis when it propagates along the line-of-sight path. The mean is 0 and the variance is less than 3 degrees. It is a random variable that obeys uniform distribution or Gaussian distribution. x (t) represents the Gaussian noise vector of the uniform linear array output on the x-axis at the tth moment, n z (t) represents the Gaussian noise vector output by the uniform linear array on the z-axis at the tth moment, n x (t) and n z (t) are not related to each other.

[0055] Step 3: Using statistical averaging, perform cross-correlation operations on the N-time sampling data of the uniform linear array output on the x-axis and the N-time sampling data of the first array element output on the z-axis to obtain the cross-covariance vector r xz Similarly, using statistical averaging, the cross-correlation operation is performed on the N-time sampling data of the uniform linear array output on the z-axis and the N-time sampling data of the first array element output on the x-axis to obtain the cross-covariance vector r zx .

[0056] In this embodiment, Among them, z 1 (t) represents the sampling data of the first array element output on the z-axis at the tth time, x 1 (t) represents the sampled data of the data output by the first array element on the x-axis at the tth moment, the superscript “H” represents the conjugate transpose operation, and p represents the power vector of the received signal on the L-type array after the signals emitted by all incoherent distribution sources at N moments propagate through all their respective transmission paths, and p = [p 1 ,…,p K ] T , p k is the kth element of p, p k The signal emitted by the kth incoherent distribution source at N moments is transmitted by L k The power of the received signal on the L-shaped array after propagating through the transmission path, The superscript “*” indicates a conjugation operation.

[0057] Step 4: Using the conjugate symmetry strategy, construct r xz The cross-covariance vector of the aperture expansion array is r e ; Similarly, using the conjugate symmetry strategy, construct r zx The cross-covariance vector of the aperture expansion array is r f .

[0058] In this embodiment, Among them, r e and r fThe dimensions of are (2M-1)×1, J represents a permutation matrix in which the elements on the anti-diagonal are 1 and the other elements are 0. Indicated by r xz The subvector consisting of the 2nd to Mth elements of Indicated by r zx The subvector consisting of the 2nd to Mth elements of .

[0059] Step 5: If the number of signal sources K is known a priori, then according to r e and r f , and use a closed-form solution algorithm based on rotation-invariant subspace technology to obtain K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis.

[0060] Further limit, according to r e and r f , and the specific process of obtaining K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis using a closed-form solution algorithm based on the rotation-invariant subspace technology is as follows:

[0061] Step 5.1a: Set r e Divide into M sub-vectors with overlapping elements, and transform r e The mth sub-vector obtained by partitioning Expressed as Similarly, r f Divide into M sub-vectors with overlapping elements, and transform r f The mth sub-vector obtained by partitioning Where m = 1, ..., M, Indicated by r e The subvector consisting of the mth element to the M+m-1th element of Indicated by r f The subvector consisting of the m-th element to the M+m-1-th element of .

[0062] Here, according to Knowable Indicated by r e The subvector consisting of the 1st element to the Mth element of Indicated by r e The subvector consisting of the 2nd element to the M+1th element, and so on, Indicated by r e The subvectors from the Mth element to the 2M-1th element of each adjacent subvector have M-1 identical elements; according to Knowable Indicated by r f The subvector consisting of the 1st element to the Mth element of Indicated by r f The subvector consisting of the 2nd element to the M+1th element, and so on, Indicated by r f The subvector consisting of the Mth element to the 2M-1th element of , there are M-1 identical elements in every two adjacent subvectors.

[0063] Step 5.2a: According to r e All sub-vectors obtained by partitioning are constructed r xz The associated M×M matrix Similarly, according to r f All sub-vectors obtained by partitioning are constructed r zx The associated M×M matrix Among them, Λ c0 represents the received signal power matrix with a diagonal structure, Λ c0 The dimension is K×K, Λ c0 The element in the kth row and kth column is p k , Λ c0 The rest of the elements in are 0.

[0064] Step 5.3a: Perform eigenvalue decomposition and get The eigenvectors corresponding to the K largest eigenvalues ​​of , the K eigenvectors constitute the matrix U e , then U e Divide into two sub-matrices and Likewise, Perform eigenvalue decomposition and get The eigenvectors corresponding to the K largest eigenvalues ​​of , the K eigenvectors constitute the matrix U f , then U f Divide into two sub-matrices and Among them, U e and U f The dimensions of are all M×K, Indicated by U e The submatrix consisting of the 1st to M-1th rows of Indicated by U e The submatrix composed of the 2nd to Mth rows of Indicated by U f The submatrix consisting of the 1st to M-1th rows of Indicated by U f The submatrix consisting of the 2nd to Mth rows of .

[0065] Step 5.4a: Calculation Then for He Perform eigenvalue decomposition to obtain K eigenvalues ​​u 1 , ..., u K ; Then through the formula k'∈[1, K] calculates θ k' , the estimated value of k'∈[1,K] k'∈[1,K]; where superscript represents the pseudo-inverse operation of the matrix, u k' Expressing e The k'th eigenvalue obtained by eigenvalue decomposition, ∠u k' Indicates u k' Phase operation, θ k' represents the k'th direction angle based on the x-axis, Represents θ k' The estimated value, however It is not necessarily the direction angle θ of the signal emitted by the kth incoherent distribution source based on the x-axis. k The estimated value of .

[0066] Similarly, calculation Then for H f Perform eigenvalue decomposition to obtain K eigenvalues ​​v 1 , ..., v K , and then through the formula k'∈[1, K] calculates β k' , the estimated value of k'∈[1,K] k'∈[1,K]; where v k' Expressing f The k'th eigenvalue obtained by eigenvalue decomposition, ∠v k' Indicates v k' Phase operation, β k' represents the k'th direction angle based on the z-axis, Represents β k' The estimated value, however It is not necessarily the direction angle β of the signal emitted by the kth incoherent distribution source based on the z axis. k The estimated value of .

[0067] If the number of signal sources K is unknown, then according to r e and r f , and use the sparse total least squares algorithm to obtain K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis.

[0068] Further limit, according to r e and r f, and the specific process of using the sparse total least squares algorithm to obtain K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis is:

[0069] Step 5.1b: Set r e Sparse representation is r e =(Ψ(θ)+E θ ) Nθ +e θ , r f Sparse representation is r f =(Ψ(β)+E β ) Nβ +e β ; where Ψ(θ) represents the matrix [JA * (θ), A (2:M) (θ)] T The sparse representation of A (2:M) (θ) represents the submatrix consisting of the 2nd to Mth rows of A(θ), Ψ(β) represents the matrix [JA * (β), A (2:M) (β)] T The sparse representation of A (2:M) (β) represents the submatrix consisting of the 2nd to Mth rows of A(β), E θ represents the θ-related deviation matrix caused by the application of generalized array flow pattern, E β represents the β-related deviation matrix caused by the application of generalized array flow pattern, p Nθ and p Nβ are all K sparse vectors, p Nθ and p Nβ contains both non-zero elements and zero elements, p Nθ The index value of the k'th non-zero element in is θ 1 , ..., θ K One of the Nβ The index value of the k'th non-zero element in is β 1 , ..., β K One of them, k'∈[1,K], that is, p Nθ The positions of the non-zero elements in represent the 1 , ..., θ K Non-one-to-one corresponding arrival angle information, p Nβ The position of the non-zero elements in represents the 1 , ..., β K Non-one-to-one corresponding arrival angle information, e θ Represents r e The estimated value and the true value (that is, the sampling data at N moments tend to infinity) e The deviation vector between the estimated values ​​of βRepresents r f The estimated value and the true value (that is, the sampling data at N moments tend to infinity) f The deviation vector between the estimated values ​​of θ and e β are all deviation vectors caused by finite sample calculations.

[0070] Step 5.2b: Construct the sparse total least squares optimization problem P1 associated with θ, described as: P1: Similarly, construct the β-related sparse total least squares optimization problem P2, which is described as: P2: Among them, |||| 2 It means to find the l of the vector 2 norm,|||| F It means to find the Frobenius norm of the matrix, |||| 1 It means to find the l of the vector 1 norm, λ represents the regularization parameter, λ>0, in this embodiment, λ=2, Indicates p Nθ The estimated value of Indicates E θ The estimated value of Indicates p Nβ The estimated value of Indicates E β The estimated value of .

[0071] Step 5.3b: Based on the θ-related sparse total least squares optimization problem P1 and the β-related sparse total least squares optimization problem P2, use the alternating descent suboptimal algorithm to complete p through alternating iterations Nθ and p Nβ The specific process is as follows:

[0072] Step 5.3.1b: Let iter represent the number of iterations, and the initial value of iter is 1; let iter max Indicates the maximum number of iterations. In this embodiment, iter max =3.

[0073] Step 5.3.2b: Given and Calculate E at the iterth iteration θ Estimated value of and E β Estimated value of Among them, when iter=1 For p Nθ The initial estimate of For p Nβ The initial estimate of When iter>1 Indicates that p is the iter-1th iteration Nθ The estimated value of Indicates that p is the iter-1th iteration Nβ The estimated value of .

[0074] Here, The solution process is: construct l 2 -l 1 The norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get The solution process is: construct l 2 -l 1 The norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get

[0075] Step 5.3.3b: Given and Calculate p at the iterth iteration Nθ Estimated value of and p Nβ Estimated value of The calculation process is: construct l 2 -l 1 The norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get The calculation process is: construct l 2 -l 1 The norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get

[0076] Step 5.3.4b: Determine whether iter is equal to iter max If yes, the alternating iteration process ends and the obtained and Corresponding to p Nθ The final estimate of p Nβ final estimated value; otherwise, let iter = iter + 1, and return to step 5.3.2b to continue execution; where the "=" in iter = iter + 1 is the assignment symbol.

[0077] Step 5.4b: By finding p Nθ The indices of the non-zero elements in the final estimate of θ determine1 , ..., θ K Estimated value of p Nθ The index value of the k'th non-zero element in the final estimate of is Represents the k'th direction angle θ based on the x-axis k' The estimated value, however It is not necessarily the direction angle θ of the signal emitted by the kth incoherent distribution source based on the x-axis. k Similarly, by finding p Nβ The indices of the nonzero elements in the final estimate of β determine 1 , ..., β K Estimated value of p Nβ The index value of the k'th non-zero element in the final estimate of is Indicates the k'th direction angle β based on the z axis k' The estimated value, however It is not necessarily the direction angle β of the signal emitted by the kth incoherent distribution source based on the z axis. k The estimated value of .

[0078] Step 6: Pair the K direction angle estimates based on the x-axis and the K direction angle estimates based on the z-axis to obtain K pairs of two-dimensional arrival angle estimates. If K = 1, no pairing is required, but in actual application scenarios, K is generally greater than 1.

[0079] Further defined, the specific process of step 6 is:

[0080] Step 6.1: Calculate the L-array cross-covariance matrix R xz , Among them, Λ c0 represents the received signal power matrix with a diagonal structure, Λ c0 The dimension is K×K, Λ c0 The element in the kth row and kth column is p k .

[0081] Step 6.2: R xz Perform column vectorization to obtain vector r v .

[0082] Step 6.3: According to r v , use the greedy algorithm to unambiguously pair the K direction angle estimates based on the x-axis and the K direction angle estimates based on the z-axis, as follows:

[0083] Step 6.3.1: Construct the complete basis matrix Ξ,

[0084]

[0085] ;in, represents the Kronecker product, Replace a(θ k ) k get Replace a(β k ) k get

[0086] Step 6.3.2: Calculate |Ξ H r v |, get K K × 1 dimensional output results Among them, || means taking the absolute value operation, Represents the output of the k'th K×1 dimension.

[0087] Step 6.3.3: The maximum value of k' and Complete the pairing process, if Maximum, then and is a pair of two-dimensional arrival angle estimates for the same incoherently distributed source.

[0088] The performance of the method of the present invention is analyzed through simulation experiments below.

[0089] The simulation was performed using MATLAB software. The two-dimensional arrival angles of the two incoherent distribution sources are {θ 1 =79.6°, β 1 =98.8°}, {θ 2 =105.3°, β 2 =127.2°}, the signals emitted by the two incoherent distribution sources are along L 1 =L 2 = 100 transmission paths, complex path gain γ 1，l (t), γ 2，l (t) all satisfy the Gaussian distribution with mean 0 and variance 1. The L-type array is a uniform linear array consisting of 10 array elements on the x-axis and z-axis with a spacing of half the carrier wavelength. The number of sampling snapshots in the simulation is 200. Except for simulation experiment 3, the variance of the random variables generated by the angle expansion is 1.5°.

[0090] Simulation experiment 1: This experiment is used to demonstrate the effect of the two-dimensional arrival angle estimation method of the present invention under the condition of known number of signal sources and uniform white noise, and is measured by the root mean square error index of the two-dimensional arrival angle estimation. The signal-to-noise ratio changes from -9 dB to 6 dB, and the simulation results are as follows Figure 3As shown, it can be found that the method of the present invention can obtain significantly improved estimation performance compared with the rotationally invariant subspace algorithm proposed by S. Shahbazpanahi et al. of the University of Ontario Institute of Technology in Canada, and also obtains certain performance improvements compared with the reduced rank subspace algorithm proposed by F. Gao et al. of Tsinghua University; at the same time, since the method of the present invention is a closed-form solution algorithm under the condition that the number of signal sources is known, the computational complexity is much lower than that of the reduced rank subspace algorithm.

[0091] Simulation experiment 2: This experiment is used to demonstrate the effect of the two-dimensional arrival angle estimation method of the present invention under the condition of unknown number of signal sources and non-uniform noise. The signal-to-noise ratio changes from -9 dB to 6 dB. The covariance matrix of the Gaussian distribution noise vector output by the uniform linear array on the x-axis is Q x =diag(2.0, 8.0, 1.5, 2.5, 7.0, 1.0, 5.5, 9.0, 1.0, 2.0), the covariance matrix of the Gaussian noise vector output by the uniform linear array on the z-axis is Q y =diag(2.3, 1.5, 2.0, 3.0, 5.9, 3.5, 9.0, 1.0, 3.0, 6.8), diag() is to find the diagonal matrix, the simulation results are as follows Figure 4 As shown, it can be seen that the method of the present invention can obtain an estimation performance that is significantly better than the comparison method under the condition of a signal-to-noise ratio of <0 dB, demonstrating the superior anti-noise robustness of the method of the present invention.

[0092] Simulation experiment 3: This experiment is used to demonstrate the effect of the method of the present invention on 2D arrival angle estimation under the conditions of unknown number of signal sources, non-uniform noise and different angle spread variances. The signal-to-noise ratio is fixed at -3 dB, and the covariance matrix of the Gaussian distribution noise vector output by the uniform linear array on the x-axis is Q x =diag(2.0, 8.0, 1.5, 2.5, 7.0, 1.0, 5.5, 9.0, 1.0, 2.0), the covariance matrix of the Gaussian noise vector output by the uniform linear array on the z-axis is Q y =diag(2.3, 1.5, 2.0, 3.0, 5.9, 3.5, 9.0, 1.0, 3.0, 6.8). The variance of the random variable generated by the angle expansion changes from 0.3° to 1.8°. The simulation results are as follows: Figure 5 As shown, it can be seen that as the variance of the random variable generated by the angle expansion gradually increases, the estimation performance of the method of the present invention does not decline significantly, and the estimation effect is much better than the comparison method, which proves the good applicability of the method of the present invention under the conditions of unknown number of signal sources, non-uniform noise and different angle expansion variances.

Claims

1. A robust two-dimensional arrival angle estimation method for incoherently distributed sources in an L-type array, characterized in that Here are the steps: Step 1: Establish a generalized array flow receiving model under an L-type array: assume that K far-field narrowband unrelated and incoherent distributed sources are incident on an L-type array consisting of two uniform linear arrays. The two uniform linear arrays are respectively located on the x-axis and z-axis of the xoz plane. The number of array elements is M, M>K, and the array element spacing is 1 / 2 of the carrier wavelength. The two uniform linear arrays have no shared array elements. The array element closest to the coordinate origin of the xoz plane in the uniform linear array on the x-axis is taken as the first array element on the x-axis, and the array element closest to the coordinate origin of the xoz plane in the uniform linear array on the z-axis is taken as the first array element on the z-axis. Step 2: In the generalized array flow receiving model, the data output by the uniform linear array on the x-axis is sampled at N times; similarly, the data output by the uniform linear array on the z-axis is sampled at N times; Step 3: Using statistical averaging, perform cross-correlation operations on the N-time sampling data of the uniform linear array output on the x-axis and the N-time sampling data of the first array element output on the z-axis to obtain the cross-covariance vector r xz Similarly, using statistical averaging, the cross-correlation operation is performed on the N-time sampling data of the uniform linear array output on the z-axis and the N-time sampling data of the first array element output on the x-axis to obtain the cross-covariance vector r zx ; Step 4: Using the conjugate symmetry strategy, construct r xz The cross-covariance vector of the aperture expansion array is r e ; Similarly, using the conjugate symmetry strategy, construct r zx The cross-covariance vector of the aperture expansion array is r f ; Step 5: If the number of signal sources K is known a priori, then according to r e and r f , and use a closed-form solution algorithm based on rotation-invariant subspace technology to obtain K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis; If the number of signal sources K is unknown, then according to r e and r f , and use the sparse total least squares algorithm to obtain K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis; Step 6: Pair the K direction angle estimates based on the x-axis and the K direction angle estimates based on the z-axis to obtain K pairs of two-dimensional arrival angle estimates.

2. The method for estimating two-dimensional arrival angle of a robust incoherently distributed source under an L-shaped array according to claim 1, characterized in that In step 2, the sampled data x(t) at the tth moment of the data output by the uniform linear array on the x-axis is approximately expressed as: x(t)≈A(θ)c0(t)+A′(θ)c1(t)+n x (t), the sampling data z(t) at the tth moment of the data output by the uniform linear array on the z-axis is approximately expressed as: z(t)≈A(β)c0(t)+A′(β)c2(t)+n z (t); where t = 1, ..., N, θ represents the angular variable of the signal emitted by the incoherent distribution source based on the x-axis, β represents the angular variable of the signal emitted by the incoherent distribution source based on the z-axis, A(θ) represents the steering matrix of the uniform linear array on the x-axis, and A(θ) = [a(θ1), ..., a(θ K )],a(θ k ) is the kth column vector of A(θ), a(θ k ) represents θ k The associated steering vector, θ k represents the direction angle of the signal emitted by the kth incoherent distribution source based on the x-axis, k = 1, ..., K, A(β) represents the steering matrix of the uniform linear array on the z-axis, A(β) = [a(β1), ..., a(β K )],a(β k ) is the kth column vector of A(β), a(β k ) represents β k The associated steering vector, β k represents the direction angle of the signal emitted by the kth incoherent distribution source based on the z-axis, (θ k ,β k ) represents the two-dimensional arrival angle of the kth incoherently distributed source, Represents the value of a(θ k ) in the variable θ k The new vector obtained by differentiation is Indicates the value of a(β k ) in the variable β k The new vector obtained by derivation is c0(t), which represents the received signal vector of the L-type array after the signals emitted by all incoherent distribution sources at time t are propagated through all their respective transmission paths. c0(t) = [c 1,0 (t),...,c K,0 (t)] T , c k,0 (t) is the kth element in c0(t), c k,0 (t) represents the signal emitted by the kth incoherent distribution source at the tth time through L k The received signal on the L-shaped array after propagating through the transmission paths, s k (t) represents the signal emitted by the kth incoherent distribution source at the tth time, γ k,l (t) represents the complex path gain of the signal transmitted by the kth incoherent distribution source at the tth time along the lth transmission path, L k represents the number of transmission paths of the kth incoherent distribution source, the superscript "T" represents the transposition operation, c1(t) represents the received signal vector on the uniform linear array on the x-axis after the signals emitted by all incoherent distribution sources at the tth moment are propagated through all their respective transmission paths and expanded in their respective azimuths based on the x-axis, c1(t) = [c 1,1 (t),…,c K,1 (t)] T , c k,1 (t) is the kth element in c1(t), c k,1 (t) represents the signal emitted by the kth incoherent distribution source at the tth time through L k transmission paths and propagates at θ k The received signal on the uniform linear array on the x-axis after angular expansion in the direction, It represents the deviation between the azimuth angle based on the x-axis when the signal emitted by the kth incoherent distribution source propagates along the lth transmission path and the azimuth angle based on the x-axis when it propagates along the line-of-sight path. The mean is 0 and the variance is less than 3 degrees. It is a random variable that obeys uniform distribution or Gaussian distribution. c2(t) represents the received signal vector on the uniform linear array on the z-axis after the signals emitted by all incoherent distribution sources at the tth moment are propagated through all their respective transmission paths and expanded in their respective directional angles based on the z-axis. c2(t) = [c 1,2 (t),…,c K,2 (t)] T , c k,2 (t) is the kth element in c2(t), c k,2 (t) represents the signal emitted by the kth incoherent distribution source at the tth time through L k transmission paths and propagates in β k The received signal on the uniform linear array on the z-axis after angular expansion in the direction, It represents the deviation between the azimuth angle based on the z-axis when the signal emitted by the kth incoherent distribution source propagates along the lth transmission path and the azimuth angle based on the z-axis when it propagates along the line-of-sight path. The mean is 0 and the variance is less than 3 degrees. It is a random variable that obeys uniform distribution or Gaussian distribution. x (t) represents the Gaussian noise vector of the uniform linear array output on the x-axis at the tth moment, n z (t) represents the Gaussian noise vector output by the uniform linear array on the z-axis at the tth moment, n x (t) and n z (t) are not related to each other.

3. The method for estimating two-dimensional arrival angle of a robust incoherently distributed source under an L-shaped array according to claim 2, characterized in that In step 3, Wherein, z1(t) represents the sampled data of the data output by the first array element on the z-axis at the tth moment, x1(t) represents the sampled data of the data output by the first array element on the x-axis at the tth moment, the superscript "H" represents the conjugate transpose operation, and p represents the power vector of the received signal on the L-type array after the signals emitted by all incoherent distribution sources at N moments propagate through all their respective transmission paths, p=[p1,...,p K ] T , p k is the kth element of p, p k The signal emitted by the kth incoherent distribution source at N moments is transmitted by L k The power of the received signal on the L-shaped array after propagating through the transmission path, The superscript "*" indicates a conjugation operation.

4. The method for estimating two-dimensional arrival angle of a robust incoherently distributed source in an L-shaped array according to claim 3, characterized in that In step 4, Among them, r e and r f The dimensions of are (2M-1)×1, J represents a permutation matrix in which the elements on the anti-diagonal are 1 and the other elements are 0. Indicated by r xz The subvector consisting of the 2nd to Mth elements of Indicated by r zx The subvector consisting of the 2nd to Mth elements of .

5. The method for estimating two-dimensional arrival angle of a robust incoherently distributed source in an L-shaped array according to claim 4, characterized in that In step 5, according to r e and r f , and the specific process of obtaining K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis using a closed-form solution algorithm based on the rotation-invariant subspace technology is as follows: Step 5.1a: Set r e Divide into M sub-vectors with overlapping elements, and transform r e The mth sub-vector obtained by partitioning Expressed as Similarly, r f Divide into M sub-vectors with overlapping elements, and transform r f The mth sub-vector obtained by partitioning Where, m = 1, ..., M, Indicated by r e The subvector consisting of the mth element to the M+m-1th element of Indicated by r f The subvector consisting of the mth element to the M+m-1th element of ; Step 5.2a: According to r e All sub-vectors obtained by partitioning are constructed r xz The associated M×M matrix Similarly, according to r f All sub-vectors obtained by partitioning are constructed r zx The associated M×M matrix Among them, Λ c0 represents the received signal power matrix with a diagonal structure, Λ c0 The dimension is K×K, Λ c0 The element in the kth row and kth column is p k ; Step 5.3a: Perform eigenvalue decomposition and get The eigenvectors corresponding to the K largest eigenvalues ​​of , the K eigenvectors constitute the matrix U e , then U e Divide into two sub-matrices and Likewise, Perform eigenvalue decomposition and get The eigenvectors corresponding to the K largest eigenvalues ​​of , the K eigenvectors constitute the matrix U f , then U f Divide into two sub-matrices and Among them, U e and U f The dimensions of are all M×K, Indicated by U e The submatrix consisting of the 1st to M-1th rows of Indicated by U e The submatrix composed of the 2nd to Mth rows of Indicated by U f The submatrix consisting of the 1st to M-1th rows of Indicated by U f The submatrix composed of the 2nd to Mth rows of ; Step 5.4a: Calculation Then for H e Perform eigenvalue decomposition to obtain K eigenvalues ​​u1,...,u K ; Then through the formula k'∈[1,K] calculates θ k' , the estimated value of k'∈[1,K] Among them, the superscript represents the pseudo-inverse operation of the matrix, u k' Expressing e The k'th eigenvalue obtained by eigenvalue decomposition, ∠u k' Indicates u k' Phase operation, θ k' represents the k'th direction angle based on the x-axis, Represents θ k' An estimated value of Similarly, calculation Then for H f Perform eigenvalue decomposition to obtain K eigenvalues ​​v1,…,v K , and then through the formula k'∈[1,K] calculates β k' , the estimated value of k'∈[1,K] k'∈[1,K]; where v k' Expressing f The k'th eigenvalue obtained by eigenvalue decomposition, ∠v k' Indicates v k' Phase operation, β k' represents the k'th direction angle based on the z-axis, Represents β k' The estimated value of .

6. The method for estimating two-dimensional arrival angle of a robust incoherently distributed source in an L-shaped array according to claim 4, characterized in that In step 5, according to r e and r f , and the specific process of using the sparse total least squares algorithm to obtain K direction angle estimates based on the x-axis and K direction angle estimates based on the z-axis is: Step 5.1b: Set r e Sparse representation is r e =(Ψ(θ)+E θ ) Nθ +e θ , r f Sparse representation is r f =(Ψ(β)+E β ) Nβ +e β ; where Ψ(θ) represents the matrix [JA * (θ),A (2:M) (θ)] T The sparse representation of A (2:M) (θ) represents the submatrix consisting of the 2nd to Mth rows of A(θ), Ψ(β) represents the matrix [JA * (β),A (2:M) (β)] T The sparse representation of A (2:M) (β) represents the submatrix consisting of the 2nd to Mth rows of A(β), E θ represents the θ-related deviation matrix caused by the application of generalized array flow pattern, E β represents the β-related deviation matrix caused by the application of generalized array flow pattern, p Nθ and p Nβ are all K sparse vectors, p Nθ and p Nβ contains both non-zero elements and zero elements, p Nθ The index values ​​of the k'th non-zero elements in are θ1,…,θ K One of the Nβ The index values ​​of the k'th non-zero elements in are β1,...,β K one of them, k'∈[1,K], e θ Represents r e The deviation vector between the estimated value and the true value, e β Represents r f The deviation vector between the estimated value and the true value of Step 5.2b: Construct the sparse total least squares optimization problem P1 associated with θ, which is described as: Similarly, the β-related sparse total least squares optimization problem P2 is constructed and described as: Among them, || ||2 means finding the l2 norm of the vector, || || F represents the Frobenius norm of the matrix, || ||1 represents the l1 norm of the vector, λ represents the regularization parameter, λ>0, Indicates p Nθ The estimated value of Indicates E θ The estimated value of Indicates p Nβ The estimated value of Indicates E β An estimated value of Step 5.3b: Based on the θ-related sparse total least squares optimization problem P1 and the β-related sparse total least squares optimization problem P2, use the alternating descent suboptimal algorithm to complete p through alternating iterations Nθ and p Nβ The specific process is as follows: Step 5.3.1b: Let iter represent the number of iterations, and the initial value of iter is 1; let iter max Indicates the maximum number of iterations; Step 5.3.2b: Given and Calculate E at the iterth iteration θ Estimated value of and E β Estimated value of Among them, when iter=1 For p Nθ The initial estimate of For p Nβ The initial estimate of When iter>1 Indicates that p is the iter-1th iteration Nθ The estimated value of Indicates that p is the iter-1th iteration Nβ An estimated value of Step 5.3.3b: Given and Calculate p at the iterth iteration Nθ Estimated value of and p Nβ Estimated value of The calculation process is: construct The l2-l1 norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get The calculation process is: construct The l2-l1 norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get Step 5.3.4b: Determine whether iter is equal to iter max If yes, the alternating iteration process ends and the obtained and Corresponding to p Nθ The final estimate of p Nβ The final estimated value; otherwise, let iter = iter + 1, and then return to step 5.3.2b to continue execution; where the "=" in iter = iter + 1 is the assignment symbol; Step 5.4b: By finding p Nθ The indices of the non-zero elements in the final estimate of θ1,…,θ K Estimated value of p Nθ The index value of the k'th non-zero element in the final estimate of is Represents the k'th direction angle θ based on the x-axis k' Similarly, by finding p Nβ The indices of the nonzero elements in the final estimate of β1,…,β K Estimated value of p Nβ The index value of the k'th non-zero element in the final estimate of is Indicates the k'th direction angle β based on the z axis k' The estimated value of .

7. The method for estimating two-dimensional arrival angle of a robust incoherently distributed source in an L-shaped array according to claim 6, characterized in that In step 5.3.2b, The solution process is: construct The l2-l1 norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get The solution process is: construct The l2-l1 norm minimization optimization problem is described as: Then, we can solve it by second-order cone programming to get 8. The method for estimating two-dimensional arrival angle of a robust incoherently distributed source in an L-shaped array according to claim 5 or 6, characterized in that The specific process of step 6 is as follows: Step 6.1: Calculate the L-array cross-covariance matrix R xz , Among them, Λ c0 represents the received signal power matrix with a diagonal structure, Λ c0 The dimension is K×K, Λ c0 The element in the kth row and kth column is p k ; Step 6.2: R xz Perform column vectorization to obtain vector r v ; Step 6.3: According to r v , use the greedy algorithm to unambiguously pair the K direction angle estimates based on the x-axis and the K direction angle estimates based on the z-axis, as follows: Step 6.3.1: Construct the complete basis matrix Ξ, in, represents the Kronecker product; Step 6.3.2: Calculate |Ξ H r v |, get K K × 1 dimensional output results Among them, || means taking the absolute value operation, Represents the output of the k'th K×1 dimension; Step 6.3.3: The maximum value of k' and Complete the pairing process, if Maximum, then and is a pair of two-dimensional arrival angle estimates for the same incoherently distributed source.

Images