A DOA estimation method and device without prior information source number

Abstract

The present application provides a DOA estimation method and device without prior information source number, including obtaining an observation data matrix of an array, establishing an objective function with the goal of minimizing the sum of the error between the observation data matrix and the signal model matrix and the regularization term used to control the signal smoothness, iteratively solving the objective function, obtaining a final denoised signal vector and a final steering vector, obtaining a weight matrix through the final denoised signal vector, and obtaining a DOA estimate of the signal source to be estimated using the weight matrix and the final steering vector. The present invention does not require the number of information sources to be known in advance, has higher adaptability and practicality, and can also achieve higher DOA estimation accuracy. Compared with the compressed sensing method, the DOA of the source can be optimized in a continuous parameter space, does not require artificially preset angle estimation grids, avoids grid mismatch problems, has better super-resolution performance, and can distinguish two targets that are closer.

Description

Technical Field

[0001] One or more embodiments of the present specification relate to the technical field of array direction of arrival estimation, and in particular, to a DOA estimation method and device that does not require a priori information source numbers. Background Art

[0002] Super-resolution direction of arrival (DOA) estimation is a key issue in the field of array signal processing. As a traditional super-resolution method, the subspace-based DOA estimation method requires the known number of signal sources in practical applications. DOA estimation methods based on compressed sensing technology include grid-based DOA estimation methods and gridless DOA estimation methods. Grid-based DOA estimation methods such as SPICE usually rely on prior information such as the number of signal sources. Gridless DOA estimation methods such as the atomic norm minimization ANM method usually require Vandermonde decomposition to extract DOA information after the optimization is completed. The processing process relies on subspace-based methods and still requires information on the number of signal sources. Therefore, when the number of signal sources is unknown, how to estimate the direction of arrival of the signal source based on the information received by the array element is a technical problem that needs to be solved urgently. Summary of the invention

[0003] The present application describes a DOA estimation method and device that does not require a priori information source numbers, which can solve the above technical problems.

[0004] According to a first aspect, a DOA estimation method without prior information source number is provided. The method comprises:

[0005] Acquire an observation data matrix of the sensor array, wherein the observation data matrix is ​​a collection of snapshot data of a transmitted signal or a reflected signal of a signal source to be estimated received by the sensor array;

[0006] Establishing an objective function with the goal of minimizing the sum of the error between the observation data matrix and the signal model matrix and the regularization term for controlling signal smoothness, wherein the regularization term is constructed by a denoised signal vector, wherein the denoised signal vector is obtained by sparsely reconstructing the signal received by the sensor array, and the signal model matrix is ​​a denoised echo signal matrix constructed by the steering vector of the signal source to be estimated relative to the sensor array;

[0007] Iteratively solving the objective function to obtain a final denoised signal vector and a final steering vector;

[0008] The weight matrix is ​​obtained through the final denoised signal vector, and the weight matrix is ​​used to approach the inverse matrix of the covariance matrix of the signal source to be estimated. The DOA estimation value of the signal source to be estimated is obtained by using the weight matrix and the final steering vector.

[0009] In some embodiments, the objective function of the mth iteration process is:

[0010]

[0011] Where Z is the observation data matrix, Y is the signal model matrix, ||·|| F represents the Frobenius norm,

[0012] τ(tr[W m T(u m )]+tr(X)) is the regularization term, W m is the weight matrix of the mth iteration process, which is constructed by the denoised signal vector obtained in the m-1th iteration process. m ) is the first covariance matrix constructed by the denoised signal vector obtained in the m-1th iteration process, tr is the trace operation of the matrix, τ is the regularization parameter used to balance sparsity and fitting error, and X is the second covariance matrix;

[0013] is a constraint condition, used to indicate that the block matrix formed by the first covariance matrix and the second covariance matrix is ​​semi-positive definite.

[0014] In some embodiments, the iteratively solving the objective function to obtain a final denoised signal vector and a final steering vector specifically includes:

[0015] Using the denoised signal vector obtained in the previous iteration process, the weight matrix of the current iteration process is obtained;

[0016] Updating the denoised signal vector obtained in the previous iteration process, the steering vector obtained in the previous iteration process and the weight matrix of the current iteration process into the equivalent function of the objective function, wherein the equivalent function of the objective function is the RANM sparse optimization equivalent function of the objective function;

[0017] Solving the equivalent function of the updated objective function to obtain a denoised signal vector obtained in the current iteration process and a steering vector obtained in the current iteration process;

[0018] The above steps are repeated until the equivalent function of the objective function converges to the termination iteration threshold, thereby obtaining the final denoised signal vector and the final steering vector.

[0019] In some embodiments, solving the equivalent function of the updated objective function to obtain the denoised signal vector obtained in the current iteration process and the steering vector obtained in the current iteration process specifically includes:

[0020] Using the obstacle parameter and the obstacle function, the equivalent function of the objective function is converted into a KKT condition;

[0021] Solving the KKT condition, and using the solution of the KKT condition as an optimization function for the denoised signal vector;

[0022] The optimization function of the denoised signal vector is iteratively solved by using a double-loop L-BFGS to obtain the denoised signal vector obtained in the current iterative process and the steering vector obtained in the current iterative process.

[0023] In some more specific embodiments, the method of obtaining the DOA estimation value of the signal source to be estimated by using the weight matrix and the final steering vector specifically includes:

[0024] Using the weight matrix and the final steering vector a(θ), construct the Capon space spectrum Among them, a H (θ) is the final steering vector a(θ) transposed and conjugated vector;

[0025] The DOA estimation value of the signal source to be estimated is obtained through the constructed Capon spatial spectrum.

[0026] In some embodiments, the signal model matrix a(θ k ) is the steering vector,

[0027] λ is the wavelength of the signal when it propagates in space, d is the array element spacing, N is the number of array elements, K, θ k 、s k are all optimization variables, K represents the number of signal sources to be estimated, θ k Represents the angle value of the signal source to be estimated, s k Represents the source strength of the signal to be estimated.

[0028] In some more specific embodiments, the denoised signal vector obtained in the previous iteration process is used to obtain the weight matrix of the current iteration process, which specifically includes:

[0029] The denoised signal vector u obtained by the m-1th iteration process m-1 Construct the Hermitian Toeplitz matrix and obtain the first covariance matrix T(u m-1 );

[0030] Using the first covariance matrix T(u m-1 ) and the noise matrix I, the weight matrix W of the mth iteration process is obtained m , W m =(T(u m-1 )+∈I)-1 .

[0031] According to a second aspect, a DOA estimation device without prior information source number is provided. The device comprises:

[0032] A first processing module is used to obtain an observation data matrix of the sensor array, where the observation data matrix is ​​a collection of snapshot data of a transmitted signal or a reflected signal of a signal source to be estimated received by the sensor array;

[0033] a second processing module, for establishing an objective function with the goal of minimizing the sum of the error between the observation data matrix and the signal model matrix and a regularization term for controlling signal smoothness, wherein the regularization term is constructed by a denoised signal vector, the denoised signal vector is obtained by sparsely reconstructing the signal received by the sensor array, and the signal model matrix is ​​a denoised echo signal matrix constructed by a steering vector of the signal source to be estimated relative to the sensor array;

[0034] A third processing module, used for iteratively solving the objective function to obtain a final denoised signal vector and a final steering vector;

[0035] The fourth processing module is used to obtain a weight matrix through the final denoised signal vector, and the weight matrix is ​​used to approach the inverse matrix of the covariance matrix of the signal source to be estimated, and the DOA estimation value of the signal source to be estimated is obtained by using the weight matrix and the final steering vector.

[0036] According to a third aspect, a computer storage medium is provided, on which a computer program is stored. When the computer program is executed by one or more processors, the DOA estimation method without requiring a priori information source number as described in any of the above embodiments is implemented.

[0037] According to a fourth aspect, an electronic device is provided, comprising a memory and one or more processors, wherein a computer program is stored on the memory, and when the computer program is executed by the one or more processors, a DOA estimation method without requiring a priori information source number as described in any one of the above embodiments is implemented.

[0038] In the above method and device provided in the embodiments of this specification, a higher DOA estimation accuracy can be achieved, and there is no need to preset the number of sources, which has higher adaptability and practicality. Compared with the compressed sensing method, the DOA of the source can be optimized in the continuous parameter space, and there is no need to manually preset the angle estimation grid, avoiding the grid mismatch problem, having better super-resolution performance, and being able to distinguish two targets that are closer to each other. BRIEF DESCRIPTION OF THE DRAWINGS

[0039] In order to more clearly illustrate the technical solutions of the embodiments of the present application, the following briefly introduces the drawings required for use in the description of the embodiments. Obviously, the drawings described below are only some embodiments of the present application. For ordinary technicians in this field, other drawings can be obtained based on these drawings without creative work.

[0040] Figure 1 A schematic flow chart of a DOA estimation method without prior information source number provided in an embodiment of this specification is shown;

[0041] Figure 2 A schematic diagram of an array arrangement provided in an embodiment of this specification is shown;

[0042] Figure 3 A schematic flow chart of a DOA estimation method without prior information source number provided in an embodiment of this specification is shown;

[0043] Figure 4 A module schematic diagram of a DOA estimation device without prior information source number provided in an embodiment of this specification is shown;

[0044] Figure 5 A schematic diagram showing the comparison results of spatial spectra reconstructed by different algorithms provided in the embodiments of this specification;

[0045] Figure 6 A schematic diagram showing the variation curves of RMSE versus SNR for different algorithms provided in the embodiments of this specification;

[0046] Figure 7 A schematic diagram showing a comparison of the calculation time of different algorithms provided in the embodiments of this specification. DETAILED DESCRIPTION

[0047] The solution provided in this specification is described below in conjunction with the accompanying drawings.

[0048] In order to make the purpose, technical solutions and advantages of the embodiments of the present application clearer, the technical solutions in the embodiments of the present application will be described below in conjunction with the accompanying drawings.

[0049] In the description of the embodiments of the present application, words such as "exemplary", "for example" or "for example" are used to indicate examples, illustrations or descriptions. Any embodiment or design described as "exemplary", "for example" or "for example" in the embodiments of the present application should not be interpreted as being more preferred or more advantageous than other embodiments or designs. Specifically, the use of words such as "exemplary", "for example" or "for example" is intended to present related concepts in a concrete way.

[0050] In the description of the embodiments of the present application, the term "and / or" is merely a description of the association relationship of associated objects, indicating that three relationships may exist. For example, A and / or B may represent: A exists alone, B exists alone, and A and B exist at the same time. In addition, unless otherwise specified, the term "plurality" means two or more.

[0051] In addition, the terms "first" and "second" are used for descriptive purposes only and should not be understood as indicating or implying relative importance or implicitly indicating the indicated technical features. Therefore, the features defined as "first" and "second" may explicitly or implicitly include one or more of the features. The terms "include", "comprises", "has" and their variations all mean "including but not limited to", unless otherwise specifically emphasized.

[0052] It should be noted that when using traditional super-resolution methods such as subspace-based DOA estimation methods for calculation, it is necessary to meet the conditions of known number of signal sources, sufficient snapshots and high signal-to-noise ratio, and a sufficient number of snapshots are required to estimate the covariance matrix of the received signal to ensure the effectiveness of the algorithm. This leads to a sharp drop in the performance of traditional algorithms when a small number of snapshots are received.

[0053] In compressed sensing technology, grid-based DOA estimation methods such as SPICE usually require a predefined grid. If the actual DOA is not on this grid, the so-called "grid mismatch" phenomenon will occur, which will lead to deviations in the estimation results. These methods usually rely on prior information such as the number of signal sources.

[0054] The gridless atomic norm minimization (ANM) method is actually a convex relaxation of the atomic norm, which aims to solve the NP-hard problem of atomic norm minimization. The ANM method can realize sparse representation of signals in the continuous domain and can theoretically achieve infinite precision DOA estimation, especially when dealing with targets that are far away. However, when the distance between targets is too close to meet certain separation conditions, the ANM method may find it difficult to accurately distinguish them. After the optimization, Vandermonde decomposition is usually required to extract DOA information.

[0055] Currently effective methods still rely on some subspace methods, such as Root-MUSIC, which means that the ANM method still relies on the quantitative information of the signal source to some extent.

[0056] The present application only requires the snapshot data of the transmitted signal or reflected signal of the signal source to be estimated received by the sensor array, and no prior knowledge of the number of signal sources is required to obtain the DOA estimation value of the signal source.

[0057] like Figure 1As shown, in this application, the received snapshot data are arranged according to the time sequence of the incoming waves to obtain the observation data matrix Z, and the weight matrix W required for the current iteration is calculated. i , construct a RANM sparse optimization problem, use the fast interior point FIPM algorithm to solve the i-th constructed RANM sparse optimization problem, repeat the above steps until the sparsity converges to the termination threshold, and use the weighted matrix obtained in the last iteration to reconstruct the Capon space spectrum, which realizes the accurate reflection of the DOA distribution of space targets without the need for prior knowledge of the number of signal sources.

[0058] Specifically, when a far-field narrowband signal in space is incident on a uniform linear array with N array elements, the observation data matrix of L snapshot data received by the array system is represented as Z;

[0059] According to the observation data matrix Z, the DOA estimation problem is constructed:

[0060]

[0061] It should be noted that in the DOA estimation problem, the observation data matrix Z is a known parameter, and the signal model matrix Y is an optimization variable, which is a denoised echo signal matrix constructed by the steering vector of the signal source to be estimated. a(θ k ) is the steering vector of the signal source to be estimated,

[0062] λ is the wavelength of the signal when it propagates in space, d is the array element spacing, N is the number of array elements, K, θ k 、s k are all optimization variables, K represents the number of signal sources to be estimated, θ k Represents the angle value of the signal source to be estimated, s k Represents the source strength of the signal to be estimated.

[0063] Using the atomic l0 norm, the DOA estimation problem is expressed as a gridless sparse parameter optimization model:

[0064]

[0065] Assuming that the number of information sources in space is sparse, the gridless sparse parameter optimization model is equivalent to the DOA estimation problem.

[0066] The gridless sparse parameter optimization model is an atomic norm minimization problem, which is a type of non-deterministic polynomial time problem NP-hard and needs to be converted into a convex optimization problem that can be solved efficiently.

[0067] This application is based on the idea of ​​repeated iteration to enhance sparsity. In each iteration, the optimization result obtained in the previous iteration is used to generate an optimization problem similar to the ANM problem. The optimization result of each iteration is solved so that it is closer to the optimal solution of the gridless sparse parameter optimization model than the previous optimization result. By repeating weighted iteration, the performance of the optimization model is close to the performance of the atomic l0 norm, overcoming the limitation of the frequency parameter interval and improving the sparsity of the optimized Toeplitz weighted matrix.

[0068] Based on the DOA estimation problem, the optimization problem of the mth iteration process is constructed:

[0069]

[0070] Weight matrix W = (T (u m-1 )+∈I) -1 , where ∈ is a hyperparameter that is considered to be set. It is generally set to 0.5, I is the unit matrix, ∈I is a regularization term that makes the matrix inversion feasible, and on the other hand, it can be understood as a noise term in spatial spectrum recovery.

[0071] Denoised signal vector u m The vector used to represent the sparsely reconstructed denoised signal is used to form the first covariance matrix T(u m ), T represents the Toeplitz operator.

[0072] Hermitian matrix and is the optimization variable, τ is the regularization parameter used to balance sparsity and fitting error, M is the total number of iterations and N is the number of array elements.

[0073] During the first iteration, the weight matrix W = I, and the optimization problem becomes an ANM problem.

[0074] In the mth iteration process, when the weight matrix W is obtained by using the optimization result of the m-1th iteration process, the equivalent representation of the optimization problem of the mth iteration process is the original function:

[0075]

[0076] Where μ is the variable to be optimized, is the feasible domain of optimization variables, is a positive semidefinite cone.

[0077] It is a Toeplitz operator acting in the real number field, which realizes the mapping from 2N-1-dimensional real vector to N-dimensional Toeplitz matrix.

[0078] Such as u=[u0,u1,…,u 2N-2 ], T(u) is in the form

[0079]

[0080] w T =T * (W) T * Represents the dual function of T, which represents the function mapping from the N-dimensional Toeplitz matrix to the 2N-1-dimensional real vector:

[0081]

[0082] in B is an N-dimensional Toeplitz matrix, and Represents the real part and the imaginary part respectively. Obviously, about T * and T have the property

[0083] tr(T(u) H B)=T * (B) T u,(9)

[0084] Its Lagrangian function is expressed as L(μ,λ)=f(μ)-<μ,λ>, where <μ,λ> represents the inner product of μ and λ, and λ is the dual variable.

[0085] According to the dual optimization theory, the dual function of the original function can be further derived:

[0086]

[0087] Represents the dual feasible region.

[0088] In order to solve the above primitive function and dual function, the obstacle parameter t and obstacle function F(μ) are introduced:

[0089]

[0090] Where |B| represents the determinant of matrix B. The calculation of F(μ) satisfies:

[0091]

[0092] The original constrained optimization problem with added obstacles can be expressed as an unconstrained optimization problem

[0093]

[0094] According to the form of F(μ) (11), it is obvious that as t increases, t -1The closer F(μ) is to the constraint, the closer problem (13) becomes to the original problem (6).

[0095] The further improved KKT condition of the primal dual problem is calculated as:

[0096]

[0097] Further solving can get (μ (t) ,λ (t) ) is also an optimized central path about t:

[0098]

[0099] Where Λ=(T(u (t) )) -1 Y (t) , where only the denoised signal vector u (t) is a non-explicit solution. Further, the solution in equation (15) can be regarded as a vector of denoised signal u (t) The optimization problem is:

[0100]

[0101] Among them, H(u)=-log|T(u)|, h(u)=τ(w) T +τtr[Z H (T(u+0.5τe0)) -1 Z],

[0102] In this application, the optimization problem in equation (16) is solved by using a double-loop L-BFGS with variable parameters. (t) ,λ (t) )’s solution.

[0103] The double-loop L-BFGS with variable parameters has a fast convergence speed and less memory space, and uses the difference vectors stored in the past iterations to reconstruct the search direction close to the Newton direction.

[0104] Since the solution obtained by this method is not an analytical solution, the solution of the dual variable is obtained by mapping the original variable. Therefore, in order to ensure the reliability and convergence of the solution, (μ (t) ,λ (t) ) also needs to be judged Whether satisfied.

[0105] Based on this, Theorem 1 is introduced: For a real vector If and only if formula (19) holds It holds, that is, the symmetric matrix Q is semi-positive definite.

[0106]

[0107] In the formula Represents a positive semidefinite cone The dual cone of .

[0108] From Theorem 1, we can see that if we want to determine whether the real vector q is in It is only necessary to determine whether equation (19) holds true, which can be expressed in the form of FFT, greatly reducing the amount of calculation.

[0109] For the formula (11) It is also equivalent to

[0110]

[0111] You can set Or q = ξ, verify whether equation (19) is true by FFT operation, and then judge Is it true?

[0112] According to the definition of the dual function, when If it holds, the dual function satisfy

[0113] According to the KKT condition (14) and the F(μ) property (12), the relationship between the duality gap and the obstacle parameter t is further obtained:

[0114]

[0115] Therefore, when When updating, the corresponding t i+1 It can be updated using the following formula:

[0116]

[0117] Considering multiple iterations In other cases, define the lower bound function f LB , and rewrite the update method of the dual gap:

[0118] At the i-th iteration, if but

[0119]

[0120] In summary, by continuously iteratively updating and increasing the barrier parameter t, the duality gap is reduced. purpose, and then Gradually approach the optimal solution (μ ＊ , λ ＊ ) of the optimization problem (5) generated during the m-th iteration, so as to obtain the denoised signal vector u m , T(u m ) is the denoised signal vector u m reconstructed by sparse reconstruction. u m and T(u m ) are essentially the same. Due to the teoplitz property of T(u m ), matrix operations can be represented by the u m vector to simplify the operations.

[0121] Repeat the above process until the objective function value of the original function converges to the iteration termination threshold |f(μ i ) - f(μ i-1 )| < δ.

[0122] Use the obtained final weight matrix and the steering vector a(θ) of the signal source to be estimated to construct the Capon spatial spectrum where, a H (θ) is the vector obtained by transposing and conjugating the steering vector a(θ) of the signal source to be estimated, and through the constructed Capon spatial spectrum, the DOA estimate value of the signal source is obtained.

[0123] Specifically, use the angle value θ k of the signal source to be estimated obtained by iterative solution to obtain the steering vector a(θ k ) of the signal source to be estimated:

[0124]

[0125] This application does not require prior knowledge of the number of signal sources. By using the results of sparse optimization to reconstruct the Capon spatial spectrum, effective resolution of closer targets is achieved without prior knowledge.

[0126] Specifically, use the double-loop L-BFGS to solve the denoised signal vector in equation (16) The specific algorithm is as follows:

[0127] Set the update step size σ and initialize the matrix

[0128]

[0129] Assume that T sets of difference vectors have been stored during the i-th iteration process. Among them, when i < T, the storage amount is i:

[0130]

[0131]

[0132]

[0133] Specifically, the FIPM algorithm used to solve the optimization problem in the mth iteration is as follows:

[0134]

[0135] This embodiment proposes a DOA estimation method without prior information source number, which is used to realize fast super-resolution DOA estimation without prior information source number. Based on the idea of ​​repeated iterative enhanced sparsity, in each iteration process, the optimization result obtained in the previous iteration process is used to generate an optimization problem similar to the ANM problem, and the optimization result of each iteration is solved so that the optimization result of each iteration is closer to the optimal solution of the gridless sparse parameter optimization model than the previous optimization result. By repeating weighted iteration, the performance of the optimization model is close to the performance of the atomic l0 norm, overcoming the limitation of the frequency parameter interval, and improving the sparsity of the optimized Toeplitz weighted matrix. In addition, without the need for prior knowledge of the number of information sources, the Capon space spectrum is reconstructed using the results of sparse optimization, and the effective resolution of closer targets is achieved. The performance exceeds the ANM method and is outstanding in terms of DOA estimation accuracy and resolution. In addition, in view of the high computational complexity of the RANM sparse optimization process, a fast interior point method FIPM algorithm is derived, and a fast optimization solution based on RANM sparse is given.

[0136] Figure 3 A flow chart of a DOA estimation method without prior information source number provided in an embodiment of this specification is shown, and the method comprises the following steps:

[0137] Step 110: Obtain an observation data matrix of the sensor array. The observation data matrix is ​​a collection of snapshot data of the transmitted signal or reflected signal of the signal source to be estimated received by the sensor array.

[0138] Exemplarily, a far-field narrowband signal in space is incident on a uniform linear array having N array elements, and the snapshot data received by the array elements are arranged as an observation data matrix Z.

[0139] In some embodiments, the DOA estimation problem is constructed using the observation data matrix Z:

[0140]

[0141] Among them, Y is the signal model matrix, which is used to represent the ideal noise-free echo signal matrix,

[0142] λ is the wavelength of the signal when it propagates in space, d is the array element spacing, N is the number of array elements, K, θ k 、s k are all optimization variables, K represents the number of information sources to be calculated, θ k represents the angle estimate to be calculated, s k represents the signal source strength to be calculated, ‖·‖ F represents the Frobenius norm, min represents the minimum function, st represents the constraints that must be met, σ 2 It's energy.

[0143] In step 120, an objective function is established with the goal of minimizing the sum of the error between the observation data matrix and the signal model matrix and the regularization term for controlling signal smoothness. The regularization term is constructed by a denoised signal vector, which is obtained by sparsely reconstructing the signal received by the sensor array, and the signal model matrix is ​​a denoised echo signal matrix constructed by the steering vector of the signal source to be estimated relative to the sensor array.

[0144] Specifically, the objective function of the mth iteration process is:

[0145]

[0146] Where Z is the observation data matrix, Y is the signal model matrix, ||·|| F represents the Frobenius norm,

[0147] τ(tr[W m T(u m )]+tr(X)) is the regularization term, W m is the weight matrix of the mth iteration process, which is constructed by the denoised signal vector obtained in the m-1th iteration process. m ) is the first covariance matrix constructed from the denoised signal vector obtained in the m-1th iteration process, tr is the trace operation of the matrix, τ is the regularization parameter used to balance sparsity and fitting error, and X is the second covariance matrix.

[0148] is a constraint condition, which is used to indicate that the block matrix formed by the first covariance matrix and the second covariance matrix is ​​semi-positive definite.

[0149] In step 130, the objective function is iteratively solved to obtain a final denoised signal vector and a final steering vector.

[0150] In some embodiments, step 130 specifically includes:

[0151] 121. Using the denoised signal vector obtained in the previous iteration process, the weight matrix of the current iteration process is obtained.

[0152] Specifically, the denoised signal vector u obtained by the m-1th iteration process is m-1 Construct the Hermitian Toeplitz matrix and get the first covariance matrix T(u m-1 );

[0153] Using the first covariance matrix T(u m-1 ) and the noise matrix I, the weight matrix W of the mth iteration process is obtained m , W m =(T(u m-1 )+∈I) -1 .

[0154] 122. The denoised signal vector obtained in the previous iteration process, the steering vector obtained in the previous iteration process and the weight matrix of the current iteration process are updated into the equivalent function of the objective function, where the equivalent function of the objective function is the RANM sparse optimization equivalent function of the objective function.

[0155] Among them, the equivalent functions of the objective function include the primitive function and the dual function. The primitive function is the equivalent function representation of the objective function, and the dual function is derived from the primitive function.

[0156] Specifically, in the mth iteration process, when the weight matrix W is obtained using the optimization result of the m-1th iteration process, the equivalent representation of the optimization problem of the mth iteration process is the original function:

[0157]

[0158] Where μ is the variable to be optimized, is the feasible domain of optimization variables, is a positive semidefinite cone.

[0159] It is a Toeplitz operator acting in the real number field, which realizes the mapping from 2N-1-dimensional real vector to N-dimensional Toeplitz matrix.

[0160] Such as u=[u0,u1,…,u 2N-2 ], T(u) is in the form

[0161]

[0162] w T =T * (W) T * Represents the dual function of T, which represents the function mapping from the N-dimensional Toeplitz matrix to the 2N-1-dimensional real vector:

[0163]

[0164] in B is an N-dimensional Toeplitz matrix, and Represents the real part and the imaginary part respectively. Obviously, about T * and T have the property

[0165] tr(T(u) H B)=T * (B) T u(9)

[0166] Its Lagrangian function is expressed as L(μ,λ)=f(μ)-<μ,λ>, where <μ,λ> represents the inner product of μ and λ, and λ is the dual variable.

[0167] According to the dual optimization theory, the dual function of the original function can be further derived:

[0168]

[0169] Represents the dual feasible region.

[0170] 123. Solve the equivalent function of the updated objective function to obtain the denoised signal vector obtained in the current iteration process and the steering vector obtained in the current iteration process.

[0171] In some embodiments, step 123 specifically includes:

[0172] 1241. Using obstacle parameters and obstacle functions, the equivalent function of the objective function is transformed into KKT conditions.

[0173] 1242. Solve the KKT condition, and use the solution of the KKT condition as the optimization function for the denoised signal vector.

[0174] 1243. Double-loop L-BFGS is used to iteratively solve the optimization function of the denoised signal vector to obtain the denoised signal vector obtained in the current iteration process and the steering vector obtained in the current iteration process.

[0175] 124. Repeat the above steps until the equivalent function of the objective function converges to the termination iteration threshold, and obtain the final denoised signal vector and the final steering vector.

[0176] The optimization variable θ obtained by iterative solution k Enter the following formula:

[0177] The steering vector a(θ) of the signal source to be estimated is obtained.

[0178] Step 140, obtain a weight matrix through the final denoised signal vector, the weight matrix is ​​used to approach the inverse matrix of the covariance matrix of the signal source to be estimated, and use the weight matrix and the final steering vector to obtain the DOA estimation value of the signal source to be estimated.

[0179] In some embodiments, the denoised signal vector u obtained by the m-1th iteration process is m-1 Construct the HermitianToeplitz matrix and get the first covariance matrix T(u m-1 );

[0180] Using the first covariance matrix T(u m-1 ) and the noise matrix I, the weight matrix W of the mth iteration process is obtained m , W m =(T(u m-1 )+∈I) -1 .

[0181] In some embodiments, a weight matrix is ​​used And the final steering vector a(θ), construct the Capon space spectrum Among them, a H (θ) is the vector obtained by transposing and conjugating the final steering vector a(θ) of the signal source to be estimated;

[0182] The DOA estimation value of the signal source to be estimated is obtained by constructing the Capon space spectrum.

[0183] The above embodiment can achieve higher DOA estimation accuracy, and does not require the number of sources to be preset, and has higher adaptability and practicality. Compared with the compressed sensing method, the DOA of the source can be optimized in the continuous parameter space, and there is no need to manually preset the angle estimation grid, avoiding the grid mismatch problem, having better super-resolution performance, and being able to distinguish two targets that are closer to each other.

[0184] Corresponding to the above method provided by this application, Figure 4 The structure diagram of a DOA estimation device without prior information source number provided by an embodiment of this specification is shown. The device includes:

[0185] A first processing module is used to obtain an observation data matrix of the sensor array, where the observation data matrix is ​​a collection of snapshot data of the transmitted signal or reflected signal of the signal source to be estimated received by the sensor array;

[0186] Exemplarily, a far-field narrowband signal in space is incident on a uniform linear array having N array elements, and the snapshot data received by the array elements are arranged as an observation data matrix Z.

[0187] In some embodiments, the DOA estimation problem is constructed using the observation data matrix Z:

[0188]

[0189] Among them, Y is the signal model matrix, which is used to represent the ideal noise-free echo signal matrix,

[0190] λ is the wavelength of the signal when it propagates in space, d is the array element spacing, N is the number of array elements, K, θ k 、s k are all optimization variables, K represents the number of information sources to be calculated, θ k represents the angle estimate to be calculated, s k represents the signal source strength to be calculated, ‖·‖ F represents the Frobenius norm, min represents the minimum function, st represents the constraints that must be met, σ 2 It's energy.

[0191] The second processing module is used to establish an objective function with the goal of minimizing the sum of the error between the observation data matrix and the signal model matrix and the regularization term used to control the signal smoothness.

[0192] The regularization term is constructed by a denoised signal vector, which is obtained by sparsely reconstructing the signal received by the sensor array. The signal model matrix is ​​a denoised echo signal matrix constructed by the steering vector of the signal source to be estimated relative to the sensor array.

[0193] Specifically, the objective function of the mth iteration process is:

[0194]

[0195] Where Z is the observation data matrix, Y is the signal model matrix, ||·|| F represents the Frobenius norm,

[0196] τ(tr[W m T(u m )]+tr(X)) is the regularization term, W m is the weight matrix of the mth iteration process, which is constructed by the denoised signal vector obtained in the m-1th iteration process. m ) is the first covariance matrix constructed from the denoised signal vector obtained in the m-1th iteration process, tr is the trace operation of the matrix, τ is the regularization parameter used to balance sparsity and fitting error, and X is the second covariance matrix.

[0197] is a constraint condition, which is used to indicate that the block matrix formed by the first covariance matrix and the second covariance matrix is ​​semi-positive definite.

[0198] The third processing module is used to iteratively solve the objective function to obtain the final denoised signal vector and the final steering vector.

[0199] In some embodiments, the third processing module specifically includes:

[0200] The first processing unit is used to obtain a weight matrix of a current iteration process by using a denoised signal vector obtained in a previous iteration process.

[0201] Specifically, the denoised signal vector u obtained by the m-1th iteration process is m-1 Construct the Hermitian Toeplitz matrix and get the first covariance matrix T(u m-1 );

[0202] Using the first covariance matrix T(u m-1 ) and the noise matrix I, the weight matrix W of the mth iteration process is obtained m , W m =(T(u m-1 )+∈I) -1 .

[0203] The second processing unit is used to update the denoised signal vector obtained in the previous iteration process, the steering vector obtained in the previous iteration process and the weight matrix of the current iteration process into the equivalent function of the objective function, and the equivalent function of the objective function is the RANM sparse optimization equivalent function of the objective function.

[0204] Among them, the equivalent functions of the objective function include the primitive function and the dual function. The primitive function is the equivalent function representation of the objective function, and the dual function is derived from the primitive function.

[0205] Specifically, in the mth iteration process, when the weight matrix W is obtained using the optimization result of the m-1th iteration process, the equivalent representation of the optimization problem of the mth iteration process is the original function:

[0206]

[0207] Where μ is the variable to be optimized, is the feasible domain of optimization variables, is a positive semidefinite cone.

[0208] It is a Toeplitz operator acting in the real number field, which realizes the mapping from 2N-1-dimensional real vector to N-dimensional Toeplitz matrix.

[0209] Such as u=[u0,u1,…,u 2N-2 ], T(u) is in the form

[0210]

[0211] w T =T * (W) T * Represents the dual function of T, which represents the function mapping from the N-dimensional Toeplitz matrix to the 2N-1-dimensional real vector:

[0212]

[0213] in B is an N-dimensional Toeplitz matrix, and Represents the real part and the imaginary part respectively. Obviously, about T * and T have the property

[0214] tr(T(u) H B)=T * (B) T u(9)

[0215] Its Lagrangian function is expressed as L(μ,λ)=f(μ)-<μ,λ>, where <μ,λ> represents the inner product of μ and λ, and λ is the dual variable.

[0216] According to the dual optimization theory, the dual function of the original function can be further derived:

[0217]

[0218] Represents the dual feasible region.

[0219] The third processing unit is used to solve the equivalent function of the updated objective function to obtain the denoised signal vector obtained in the current iterative process and the steering vector obtained in the current iterative process.

[0220] In some embodiments, the third processing unit specifically includes:

[0221] The first sub-processing unit is used to convert the equivalent function of the objective function into a KKT condition by using the obstacle parameter and the obstacle function.

[0222] The second sub-processing unit is used to solve the KKT condition and use the solution of the KKT condition as the optimization function about the denoised signal vector.

[0223] The third sub-processing unit is used to iteratively solve the optimization function about the denoised signal vector by using double-loop L-BFGS to obtain the denoised signal vector obtained in the current iterative process and the steering vector obtained in the current iterative process.

[0224] The fourth processing unit is used to repeat the above steps until the equivalent function of the objective function converges to the termination iteration threshold, thereby obtaining the final steering vector and the final steering vector.

[0225] The optimization variable θ obtained by iterative solution k Enter the following formula:

[0226] The steering vector a(θ) of the signal source to be estimated is obtained.

[0227] The fourth processing module is used to obtain a weight matrix through the final denoised signal vector. The weight matrix is ​​used to approach the inverse matrix of the covariance matrix of the signal source to be estimated. The DOA estimation value of the signal source to be estimated is obtained using the weight matrix and the final vector to be guided.

[0228] In some embodiments, in some embodiments, the denoised signal vector u obtained by the m-1th iteration process is m-1 Construct the Hermitian Toeplitz matrix and get the first covariance matrix T(u m-1 );

[0229] Using the first covariance matrix T(u m-1 ) and the noise matrix I, the weight matrix W of the mth iteration process is obtained m , W m =(T(u m-1 )+∈I) -1 .

[0230] In some embodiments, a weight matrix is ​​used And the final steering vector a(θ), construct the Capon space spectrum Among them, a H (θ) is the vector obtained by transposing and conjugating the final steering vector a(θ).

[0231] The DOA estimation value of the signal source to be estimated is obtained by constructing the Capon space spectrum.

[0232] According to another embodiment, there is also provided a computer-readable storage medium having a computer program stored thereon, which, when executed in a computer, causes the computer to execute a combination of Figure 3 The method described.

[0233] According to another embodiment of the present invention, a computing device is provided, including a memory and a processor, wherein the memory stores executable code, and when the processor executes the executable code, the Figure 3 The method described.

[0234] The following are experimental results obtained by using numerical simulation experiments to verify the method of the present application and analyze it.

[0235] In order to verify the performance of the algorithm proposed in this application, numerical simulation experiments are used for analysis. Consider a sonar array system with 20 array elements, the number of received snapshots L = 2, the signal-to-noise ratio (SNR) is 25dB, the targets are set at 5°, 9°, 10°, 15°, and 30°, and the signal wavelength is set to 0.15m.

[0236] like Figure 5 Figure 2 shows the performance comparison of spatial spectra drawn by several traditional algorithms and sparse algorithms. It can be seen that too small a snapshot number will lead to the inability to accurately obtain the covariance matrix, and then the MUSIC algorithm and MVDR method based on the covariance matrix to restore the spatial spectrum cannot accurately estimate the target DOA, while the sparse method can effectively deal with the problem of too few snapshots. And compared with the fast ANM (FANM) algorithm, the proposed algorithm uses sparse enhancement technology to have a stronger resolution for adjacent targets, such as Figure 5 shown.

[0237] in addition, Figure 5 The spatial spectrum performance of the RANM algorithm based on SDPT3 (RANM-SDPT3) is also demonstrated, and its resolution performance is weaker than that of the proposed algorithm in the presence of noise.

[0238] In addition, A=200 Monte Carlo experiments were performed under different signal-to-noise ratios to verify the advantages of the algorithm proposed in this application in terms of noise resistance performance.

[0239] The root mean square error (RMSE) is used as the metric, which is defined as in represents the estimated DOA of the kth target in the ath experiment, θ k Represents the true value of the kth DOA.

[0240] The RMSE estimated by the algorithm proposed in this application and the RANM-SDPT3 algorithm under different signal-to-noise ratios are compared, such as Figure 6 The results show that the algorithm proposed in this application is comparable to the RANM-SDPT3 algorithm in terms of accuracy and has good robustness and estimation accuracy.

[0241] Finally, the calculation time of RANM-SDPT3 and the algorithm proposed in this application are compared. The results are as follows: Figure 7As shown in the figure, as the number of snapshots L that the algorithm needs to process increases, the computational complexity of the optimization problem increases significantly, which creates a huge computational burden on the convex optimization tool SDPT3. However, the proposed algorithm developed based on FastIPM is less affected and has a significant advantage in computing speed compared to the RANM-SDPT3 algorithm.

[0242] Those skilled in the art should be aware that in one or more of the above examples, the functions described in this application can be implemented by hardware, software, firmware, or any combination thereof. When implemented by software, these functions can be stored in a computer-readable medium or transmitted as one or more instructions or codes on a computer-readable medium.

[0243] The specific implementation methods described above further illustrate the purpose, technical solutions and beneficial effects of the present application in detail. It should be understood that the above description is only the specific implementation methods of the present application and is not intended to limit the scope of protection of the present application. Any modifications, equivalent substitutions, improvements, etc. made on the basis of the technical solutions of the present application should be included in the scope of protection of the present application.

Claims

1. A DOA estimation method without prior information source number, characterized in that: The method comprises: Acquire an observation data matrix of the sensor array, wherein the observation data matrix is ​​a collection of snapshot data of a transmitted signal or a reflected signal of a signal source to be estimated received by the sensor array; Establishing an objective function with the goal of minimizing the sum of the error between the observation data matrix and the signal model matrix and the regularization term for controlling signal smoothness, wherein the regularization term is constructed by a denoised signal vector, the denoised signal vector is obtained by sparsely reconstructing the signal received by the sensor array, and the signal model matrix is ​​a denoised echo signal matrix constructed by the steering vector of the signal source to be estimated relative to the sensor array; Iteratively solving the objective function to obtain a final denoised signal vector and a final steering vector; The weight matrix is ​​obtained through the final denoised signal vector, and the weight matrix is ​​used to approach the inverse matrix of the covariance matrix of the signal source to be estimated. The DOA estimation value of the signal source to be estimated is obtained by using the weight matrix and the final steering vector.

2. The method according to claim 1, characterized in that: The objective function of the mth iteration process is: Where Z is the observation data matrix, Y is the signal model matrix, ||·|| F represents the Frobenius norm, τ(tr[W m T(u m )]+tr(X)) is the regularization term, W m is the weight matrix of the mth iteration process, which is constructed by the denoised signal vector obtained in the m-1th iteration process. m ) is the first covariance matrix constructed by the denoised signal vector obtained in the m-1th iteration process, tr is the trace operation of the matrix, τ is the regularization parameter used to balance sparsity and fitting error, and X is the second covariance matrix; is a constraint condition, used to indicate that the block matrix formed by the first covariance matrix and the second covariance matrix is ​​semi-positive definite.

3. The method according to claim 1, characterized in that The iteratively solving the objective function to obtain a final denoised signal vector and a final steering vector specifically includes: Using the denoised signal vector obtained in the previous iteration process, the weight matrix of the current iteration process is obtained; Updating the denoised signal vector obtained in the previous iteration process, the steering vector obtained in the previous iteration process and the weight matrix of the current iteration process into the equivalent function of the objective function, wherein the equivalent function of the objective function is the RANM sparse optimization equivalent function of the objective function; Solving the equivalent function of the updated objective function to obtain a denoised signal vector obtained in the current iteration process and a steering vector obtained in the current iteration process; The above steps are repeated until the equivalent function of the objective function converges to the termination iteration threshold, thereby obtaining the final denoised signal vector and the final steering vector.

4. The method according to claim 3, characterized in that Solving the equivalent function of the updated objective function to obtain the denoised signal vector obtained in the current iteration process and the steering vector obtained in the current iteration process specifically includes: Using the obstacle parameter and the obstacle function, the equivalent function of the objective function is converted into a KKT condition; Solving the KKT condition, and using the solution of the KKT condition as an optimization function for the denoised signal vector; The optimization function of the denoised signal vector is iteratively solved by using a double-loop L-BFGS to obtain the denoised signal vector obtained in the current iterative process and the steering vector obtained in the current iterative process.

5. The method according to claim 1, characterized in that The step of obtaining the DOA estimation value of the signal source to be estimated by using the weight matrix and the final steering vector specifically includes: Using the weight matrix and the final steering vector a(θ), construct the Capon space spectrum Among them, a H (θ) is the final steering vector a(θ) transposed and conjugated vector; The DOA estimation value of the signal source to be estimated is obtained through the constructed Capon spatial spectrum.

6. The method according to claim 1, characterized in that The signal model matrix a(θ k ) is the steering vector, λ is the wavelength of the signal when it propagates in space, d is the array element spacing, N is the number of array elements, K, θ k 、s k are all optimization variables, K represents the number of signal sources to be estimated, θ k Represents the angle value of the signal source to be estimated, s k represents the source strength of the signal source to be estimated, and T is the transpose operation of the matrix or vector.

7. The method according to claim 3, characterized in that The denoised signal vector obtained in the previous iteration process is used to obtain the weight matrix of the current iteration process, which specifically includes: The denoised signal vector u obtained by the m-1th iteration process m-1 Construct the Hermitian Toeplitz matrix and get the first covariance matrix T(u m-1 ); Using the first covariance matrix T(u m-1 ) and the noise matrix I, the weight matrix W of the mth iteration process is obtained m , W m =(T(u m-1 )+∈I) -1 .

8. A DOA estimation device without prior information source number, characterized in that: The device comprises: A first processing module is used to obtain an observation data matrix of the sensor array, where the observation data matrix is ​​a collection of snapshot data of a transmitted signal or a reflected signal of a signal source to be estimated received by the sensor array; A second processing module is used to establish an objective function with the goal of minimizing the sum of the error between the observation data matrix and the signal model matrix and the regularization term for controlling signal smoothness, wherein the regularization term is constructed by a denoised signal vector, the denoised signal vector is obtained by sparsely reconstructing the signal received by the sensor array, and the signal model matrix is ​​a denoised echo signal matrix constructed by the steering vector of the signal source to be estimated relative to the sensor array; A third processing module, used for iteratively solving the objective function to obtain a final denoised signal vector and a final steering vector; The fourth processing module is used to obtain a weight matrix through the final denoised signal vector, and the weight matrix is ​​used to approach the inverse matrix of the covariance matrix of the signal source to be estimated, and the DOA estimation value of the signal source to be estimated is obtained by using the weight matrix and the final steering vector.

9. A computer storage medium, characterized in that The computer-readable storage medium stores a computer program, and when the computer program is executed by one or more processors, the DOA estimation method without requiring a priori information source number as claimed in any one of claims 1 to 7 is implemented.

10. An electronic device, characterized in that: The method comprises a memory and one or more processors, wherein a computer program is stored in the memory, and when the computer program is executed by the one or more processors, the DOA estimation method without prior information source number as claimed in any one of claims 1 to 7 is implemented.

Images