Incoming wave direction estimation method under array sensor calibration error

Abstract

The invention discloses an incoming wave direction estimation method under array sensor calibration errors, which comprises the following steps of: 1, selecting a time domain sampling point number T according to an incident signal sampling rule, establishing an array sampling matrix system Y = B.S + E + N under the condition of sensor calibration errors, and taking the sampling point number T as a column number of a left-end item matrix, the matrix Y is an M * T array sampling matrix, M is the number of array sensors, the matrix B is an M * K array manifold matrix, K is the number of signal sources, the matrix S is a K * T signal matrix, the matrix E is an M * T calibration error matrix, and the matrix N is an M * T noise matrix; 2, an over-complete dictionary matrix A of the array manifold matrix is constructed according to a preset rule, an array sampling matrix system Y = A * X + E + N under sparse representation is obtained, the matrix A is an M * N over-complete dictionary matrix, N is the column number of a dictionary, and the matrix X is an N * T signal sparse representation matrix; 3, under a low-rank and row sparse decomposition framework, constructing a convex optimization formula according to a matrix system Y = A * X + E + N of array sampling under sparse representation; 4, solving the constructed convex optimization formula by adopting an iterative re-weighting method to obtain solution matrixes # imgabs0 # and # imgabs1 #; and 5, respectively obtaining corresponding incoming wave direction estimation and error sensor calibration positions according to non-zero line positions of the solution matrixes # imgabs2 # and # imgabs3 #. According to the method, two solution matrixes # imgabs4 # and # imgabs5 # can be calculated by adopting an iterative re-weighting method, and incoming wave direction estimation and error sensor position calibration can be respectively obtained through the solution matrixes.

Description

Technical Field

[0001] The present invention relates to the technical field of wave direction estimation, and in particular to a wave direction estimation method under calibration error of an array sensor. Background Art

[0002] In the study of wave direction estimation, compressed sensing technology is often used to solve this type of problem. The process of using compressed sensing technology to solve the wave direction estimation problem is as follows: First, the matrix form of the array sampling data can be expressed as:

[0003] Y=B·S+N

[0004] Among them, matrix Y is the M×T array sampling matrix, M is the number of array sensors, T is the number of sampling points, matrix B is the M×K array manifold matrix, K is the number of non-correlated signal sources, matrix S is the K×T signal matrix, and matrix N is the M×T noise matrix.

[0005] According to the space partitioning, the overcomplete dictionary matrix A of the array manifold matrix is ​​constructed, and the array sampling matrix is ​​rewritten as:

[0006] Y=A·X+N

[0007] Among them, the matrix A is an M×N overcomplete dictionary matrix, and the matrix X is an N×T signal sparse representation matrix, in which each non-zero row is equal to the row at the corresponding position of S. According to the array sampling matrix under sparse representation, the following convex optimization formula can be constructed to restore X:

[0008]

[0009] where ||·|| F represents the norm of matrix F, ||·|| 2，1 It means that the 2-norm of each row of the matrix is ​​calculated and then the 1-norm of the obtained vector is calculated. By solving the matrix The position of the non-zero row can be used to estimate the direction of the incoming wave.

[0010] However, the above method will not work when the array sensor is miscalibrated. The current solutions to this problem include the maximum likelihood estimation method, but this method requires that the number of correctly calibrated sensors exceeds the number of signals; rotation invariance technology is also used to jointly estimate the direction of the incoming wave, sensor gain and phase in a uniform linear array, but this method requires the location of the miscalibrated sensor to be known in advance; the maximum entropy criterion is also used to develop a wave estimation method for the unknown location of the miscalibrated sensor, but this method cannot detect the location of the miscalibrated sensor.

[0011] A solution is proposed for obtaining the direction of the incoming wave and the position of the calibrated sensor at the same time when the position of the calibrated sensor and other related information are unknown. The details are as follows. First, the matrix form of the array sampling data in the case of array sensor calibration error can be expressed as:

[0012] Y=B·S+E+N

[0013] Among them, the matrix E is the calibration error matrix of M×T. According to the space partitioning, the overcomplete dictionary matrix A of the array manifold matrix is ​​constructed, and the sparse representation of the array sampling matrix in the case of array sensor calibration error is obtained as follows:

[0014] Y=A·X+E+N

[0015] Among them, each non-zero row of the matrix X is equal to the row of the corresponding position in S, and each non-zero row of the matrix E corresponds to the position of the sensor with incorrect calibration. Based on the sparse representation of the array sampling matrix, the following convex optimization formula can be constructed:

[0016]

[0017] Where μ is the smoothing coefficient, λ1 and λ2 are the trade-off coefficients, 1 N is a column vector with N columns and all 1s, 1 M is a column vector with all 1s and M columns. The solution matrix can be obtained by iterative reweighting method and By solving the matrix and The position of the non-zero row of can be used to estimate the incoming wave direction and obtain the position of the calibrated error sensor in the same system framework. Summary of the invention

[0018] The purpose of the present invention is to provide a method for estimating the direction of arrival of an array sensor under calibration errors, by constructing an overcomplete representation of the array manifold matrix to obtain a sparse representation of the array sampling matrix; by quickly solving the row-sparse solution matrix through an iterative reweighting method; by indexing the non-zero row positions of the solution matrix, the direction of arrival of the wave is estimated and the position of the calibrated error sensor is finally obtained in the same system framework.

[0019] To achieve the above objectives, combined Figure 1 , Figure 2 The present invention proposes a method for estimating the direction of incoming waves under calibration errors of an array sensor, the method comprising:

[0020] Step 1: According to the sampling rule of the incident signal, the number of time domain sampling points T is selected, and the matrix system Y=B·S+E+N for array sampling under sensor calibration error is established. The number of sampling points T is used as the number of columns of the left-hand matrix, where the matrix Y is an M×T array sampling matrix, M is the number of array sensors, the matrix B is an M×K array manifold matrix, K is the number of signal sources, the matrix S is a K×T signal matrix, the matrix E is an M×T calibration error matrix, and the matrix N is an M×T noise matrix;

[0021] Step 2: Construct the dictionary A according to the preset rules to obtain the matrix system Y=A·X+E+N of array sampling under sparse representation, where the matrix A is an M×N overcomplete dictionary matrix, N is the number of columns of the dictionary, and the matrix X is an N×T signal sparse representation matrix;

[0022] Step 3: Under the framework of low rank and row sparse decomposition, a convex optimization formula is constructed according to the matrix system Y = A·X + E + N of array sampling under sparse representation;

[0023] Step 4: Use iterative reweighting method to solve the convex optimization formula to obtain solution matrices X and E;

[0024] Step 5: According to the non-zero row positions of the solution matrices X and E, the corresponding arrival direction estimation and the calibration error sensor position are obtained respectively.

[0025] In a further embodiment, the method further comprises:

[0026] The incident direction is set to change along the same direction, using the zenith angle θ and azimuth angle in the spherical coordinate system. represents the incident direction, where the zenith angle θ ranges from -90 to 90 degrees and the azimuth angle The range is 0-360 degrees.

[0027] In a further embodiment, in step S2, constructing the dictionary A according to a preset rule to obtain the matrix system Y=A·X+E+N under sparse representation means:

[0028] According to a certain spatial interval, the zenith angle θ and the azimuth angle Sampling is performed and an array manifold matrix over-complete dictionary A is constructed based on the sampling points.

[0029] In a further embodiment, the convex optimization formula in step S3 is as follows:

[0030]

[0031] In a further embodiment, the iterative reweighting method in step S4 comprises the following steps:

[0032] S401: Set the smoothing coefficient μ, the trade-off coefficients λ1 and λ2, and the maximum number of iterations I max and threshold ε.

[0033] S402: Initial iteration number i=0, initialize solution matrix X and E, denoted as X 0 and E 0 .

[0034] S403: Calculate the transition matrix P of the i-th iteration i and Q i , the specific calculation method is as follows:

[0035]

[0036] where ||·||2 represents the vector 2-norm, diag{·} represents the construction of a diagonal matrix using the elements in {}, and X i-1 (m) and E i-1 (m) represents the nth row of X and the mth row of E obtained by the i-1th iteration update, n = 1, 2, ..., N, m = 1, 2, ..., M.

[0037] S404: Calculate the solution matrices X and E of the i-th iteration. The specific calculation method is as follows:

[0038] X i =(A H ·A+λ1P i ) -1 ·A H ·(YE i-1 )

[0039] E i =(I+λ2Q i ) -1 ·[YA·X i ]

[0040] in(·) H represents the matrix conjugate transpose operation, (·) -1 Represents a matrix inversion operation.

[0041] S405: Calculate the relative error. The specific calculation method is as follows:

[0042]

[0043] Where f is the objective function, and the specific calculation method is as follows:

[0044]

[0045] After the calculation is completed, the number of iterations is updated i=i+1.

[0046] S406: Repeat steps S403 to S405 until i is greater than I max Or the relative error Δ is less than ε.

[0047] Compared with the prior art, the technical solution of the present invention has the following significant beneficial effects:

[0048] (1) By introducing the sensor calibration error signal into the array sampling system model, the estimation of the incoming wave direction and the acquisition of the calibration error sensor position can be performed in the same system framework.

[0049] (2) By constructing a sparse representation of the array sampling matrix in the case of array sensor calibration errors, the signal matrix and the calibration error matrix are both row-sparse, which provides a guarantee for the successful inversion of the solution matrix. In the end, the signal wave direction and the position of the array calibration error sensor are obtained without any prior information about the array and the signal source.

[0050] (3) The iterative reweighting method can be used to quickly calculate the solution matrix. According to the non-zero rows of the two solution matrices, the corresponding incoming wave direction estimation and the calibration error sensor position can be obtained respectively. BRIEF DESCRIPTION OF THE DRAWINGS

[0051] The drawings are not intended to be drawn to scale. In the drawings, each identical or nearly identical component shown in various figures may be represented by the same reference numeral. For clarity, not every component is labeled in every figure. Embodiments of various aspects of the present invention will now be described by way of example and with reference to the accompanying drawings, in which:

[0052] Figure 1 It is a flow chart of the method for estimating the direction of incoming waves under calibration error of the array sensor of the present invention;

[0053] Figure 2 is a method flow chart of the iteratively reweighted least squares method of the present invention;

[0054] Figure 3 It is the column power diagram obtained by solving the vector 2 norm of each row of the matrix X and E. DETAILED DESCRIPTION

[0055] In order to better understand the technical content of the present invention, a specific embodiment is cited and described as follows in conjunction with the attached drawings. The purpose of the present invention is to provide a method for estimating the direction of arrival of an array sensor under calibration errors, by constructing an overcomplete representation of the array manifold matrix to obtain a sparse representation of the array sampling matrix; by using an iterative reweighting method, a row-sparse solution matrix can be quickly solved; by indexing the non-zero row positions of the solution matrix, the direction of arrival of the wave and the position of the sensor with calibration errors can be estimated in the same framework.

[0056] To achieve the above objectives, combined Figure 1 , Figure 2 The present invention proposes a method for estimating the direction of incoming waves under calibration errors of an array sensor, the method comprising:

[0057] Step 1: According to the sampling rule of the incident signal, the number of time domain sampling points T is selected, and the matrix system Y=B·S+E+N for array sampling under sensor calibration error is established. The number of sampling points T is used as the number of columns of the left-hand matrix, where the matrix Y is an M×T array sampling matrix, M is the number of array sensors, the matrix B is an M×K array manifold matrix, K is the number of signal sources, the matrix S is a K×T signal matrix, the matrix E is an M×T calibration error matrix, and the matrix N is an M×T noise matrix;

[0058] Step 2: Construct the dictionary A according to the preset rules to obtain the matrix system Y=A·X+E+N under the sparse representation, where the matrix A is an M×N dictionary, N is the number of columns of the dictionary, and the matrix X is an N×T sparse signal representation;

[0059] Step 3: Under the framework of low rank and row sparse decomposition, a convex optimization formula is constructed according to the matrix system Y = A·X + E + N of array sampling under sparse representation;

[0060] Step 4: Use iterative reweighting method to solve the convex optimization formula to obtain solution matrices X and E;

[0061] Step 5: According to the non-zero row positions of the solution matrices X and E, the corresponding arrival direction estimation and the calibration error sensor position are obtained respectively.

[0062] In a specific application embodiment, the detailed steps of S1 are:

[0063] S101: Assume that the array is a linear array containing 15 sensors, and the interval between two sensors is λ / 2, where λ is the working wavelength. The signal incident direction only changes along the θ direction, and the range is [-90°, 90°]. The number of non-correlated signal sources is selected to be 4, and the incoming wave directions are θ1=-60°, θ2=-20°, θ3=20°, θ4=40°, and it is assumed that there are 50 time samples. The 5th and 10th sensors in the array are calibrated incorrectly. According to the incident signal sampling rules, the matrix system is established in the case of sensor calibration errors:

[0064] Y=B·S+E+N

[0065] Among them, matrix Y is a 15×50 array sampling matrix, matrix B is a 15×4 array manifold matrix, and each column of the matrix is ​​the steering vector corresponding to the incoming wave direction, expressed as

[0066]

[0067] where d1, d2, …, d M are the distances between the mth sensor and the first sensor, m=1, 2, ..., M, the matrix E is a 15×50 calibration error matrix, the matrix N is a 15×50 additive white Gaussian noise matrix, and the signal-to-noise ratio is 20 dB.

[0068] In a specific application embodiment, the detailed steps of S2 are:

[0069] S201: Construct a dictionary A with a 1° spatial sampling interval rule, where each column of A is a steering vector corresponding to the incoming wave direction of the sampling point, expressed as

[0070]

[0071] Thus we get the matrix system under sparse representation:

[0072] Y=A·X+E+N

[0073] Among them, matrix A is a dictionary of 15×181, and matrix X is a sparse signal representation of 181×0. Both matrix X and matrix E are row sparse matrices at this time, and this row sparse feature provides a guarantee for their solution.

[0074] In a specific application embodiment, the detailed steps of S3 are:

[0075] S301: According to the matrix system under sparse representation, the following convex optimization formula can be constructed:

[0076]

[0077] Where μ is the smoothing coefficient, λ1 and λ2 are the trade-off coefficients, 1 N is a column vector with 181 columns and all 1s. M is a column vector of all 1s with 15 columns.

[0078] In a specific application embodiment, the detailed steps of S4 are:

[0079] S401: Set smoothing coefficient μ = 10 -4 , trade-off coefficients λ1 = 8.8 and λ2 = 2, maximum number of iterations I max = 1000 and threshold ε = 10 -5 .

[0080] S402: Set the initial iteration number i=0, initialize the solution matrix X=0 and E=0, denoted as X 0 and E 0 .

[0081] S403: Calculate the transition matrix P of the i-th iteration i and Qi , the specific calculation method is as follows:

[0082]

[0083] where ||·||2 represents the vector 2-norm, diag{·} represents the construction of a diagonal matrix using the elements in {}, and X i-1 (n) and E i-1 (m) represents the nth row of X and the mth row of E obtained by the i-1th iteration update, n = 1, 2, ..., 181, m = 1, 2, ..., 15.

[0084] S404: Calculate the solution vector X of the i-th iteration i and E i , the specific calculation method is as follows:

[0085] X i =(A H ·A+λ1P i ) -1 ·A H ·(YE i-1 )

[0086] E i =(I+λ2Q i ) -1 ·[YA·X i ]

[0087] in(·) H represents the matrix conjugate transpose operation, (·) -1 Represents a matrix inversion operation.

[0088] S405: Calculate the relative error. The specific calculation method is as follows:

[0089]

[0090] Where f is the objective function, and the specific calculation method is as follows:

[0091]

[0092] After the calculation is completed, the number of iterations is updated i=i+1.

[0093] S406: Repeat steps S403 to S405 until i is greater than I max Or the relative error Δ is less than ε.

[0094] In a specific application embodiment, the detailed steps of S5 are:

[0095] S501: The column power diagram obtained by solving the vector 2 norm of each row of the matrix X and E is given by Figure 3According to the non-zero row positions of the solution matrices X and E, the corresponding estimated directions of the incoming waves are -60°, -20°, 20°, 40°, and the calibration error sensor positions are the 5th and 10th.

[0096] Various aspects of the present invention are described in this disclosure with reference to the accompanying drawings, in which many illustrative embodiments are shown. The embodiments of the present disclosure are not necessarily defined to include all aspects of the present invention. It should be understood that the various concepts and embodiments introduced above, as well as those described in more detail below, can be implemented in any of many ways, because the concepts and embodiments disclosed by the present invention are not limited to any implementation. In addition, some aspects disclosed by the present invention can be used alone, or in any appropriate combination with other aspects disclosed by the present invention.

[0097] Although the present invention has been disclosed as above with preferred embodiments, it is not intended to limit the present invention. A person with ordinary knowledge in the technical field to which the present invention belongs may make various changes and modifications without departing from the spirit and scope of the present invention. Therefore, the protection scope of the present invention shall be determined by the definition of the claims.

Claims

1. A method for estimating the direction of incoming waves when an array sensor is calibrated incorrectly, characterized in that: At the same time, the direction of the incoming wave is obtained and the position of the calibrated error sensor is corrected. The method includes: S1: According to the sampling rule of the incident signal, the number of time domain sampling points T is selected, and the matrix system Y=B·S+E+N for array sampling under sensor calibration error is established. The number of sampling points T is used as the number of columns of the left-hand matrix, where the matrix Y is an M×T array sampling matrix, M is the number of array sensors, the matrix B is an M×K array manifold matrix, K is the number of signal sources, the matrix S is a K×T signal matrix, the matrix E is an M×T calibration error matrix, and the matrix N is an M×T noise matrix; S2: construct an overcomplete dictionary matrix A of the array manifold matrix according to a preset rule, and obtain a sparse representation Y=A·X+E+N of the array sampling matrix system, where the matrix A is an M×N overcomplete dictionary matrix, N is the number of columns of the dictionary, and the matrix X is an N×T signal sparse representation matrix; S3: Construct a convex optimization formula based on the sparse representation Y = A·X+E+N of the array sampling matrix system in the framework of low rank and row sparse decomposition; S4: Use iterative reweighting to solve the convex optimization formula to obtain the solution matrix and S5: According to the solution matrix and The non-zero row positions of respectively obtain the corresponding incoming wave direction estimation and the calibration error sensor position.

2. The method for estimating the direction of arrival of a wave under an array sensor calibration error according to claim 1, characterized in that: In step S1, the sensor calibration error signal is introduced into the array sampling system model to estimate the incoming wave direction and obtain the position of the calibration error sensor.

3. The method for estimating the direction of arrival of a wave under an array sensor calibration error according to claim 1, characterized in that: In step S2, an overcomplete dictionary matrix A of the array manifold matrix is ​​constructed so that both the matrix X and the matrix E have row sparse properties.

4. A method for estimating the direction of arrival of a wave under calibration error of an array sensor according to claims 1 and 3, characterized in that: In the low-rank and row-sparse decomposition framework, the convex optimization formula described in step S3 constrains the row-sparse properties of the matrix X and the matrix E, thereby ensuring that the two solution matrices each have a unique solution.

5. The method for estimating the direction of arrival of a wave under calibration error of an array sensor according to claim 1, characterized in that: The iterative reweighting method comprises the following steps: S401: Set the smoothing coefficient μ, the trade-off coefficients λ1 and λ2, and the maximum number of iterations I max and threshold ε. S402: Initial iteration number i=0, initialize solution matrix X and E, denoted as X 0 and E 0 . S403: Calculate the transition matrix P of the i-th iteration i and Q i , the specific calculation method is as follows: where ||·||2 represents the vector 2-norm, diag{·} represents the construction of a diagonal matrix using the elements in {}, and X i-1 (n) and E i-1 (m) represents the nth row of the X matrix and the mth row of the E matrix obtained by the i-1th iteration update, n = 1, 2, ..., N, m = 1, 2, ..., M. S404: Calculate the solution matrix X of the i-th iteration i and E i , the specific calculation method is as follows: X i =(A H ·A+λ1P i ) -1 ·A H ·(Y-E i-1 ) E i =(I+λ2Q i ) -1 ·[Y-A·X i ] in(·) H represents the matrix conjugate transpose operation, (·) -1 Represents a matrix inversion operation. S405: Calculate the relative error. The specific calculation method is as follows: Where f is the objective function, and the specific calculation method is as follows: After the calculation is completed, the number of iterations is updated i=i+1. S406: Repeat steps S403 to S405 until the number of iterations i is greater than I. max Or the relative error Δ is less than ε.

Images