Two-dimensional super-resolution detection method for solving arrival direction problem of inverted halfpace array

Abstract

The invention discloses a two-dimensional super-resolution detection method for solving an arrival direction problem of an inverted halfpace array. The two-dimensional super-resolution detection method comprises the following steps: step 1, setting parameters such as a stop threshold, a maximum number of iterations, single-step and cyclic Newton correction times and the like; 2, inputting array element data Y received by the circular ring array distributed at equal intervals; 3, a detection stage: obtaining an initial estimation value of a line spectrum frequency w by using a finite discrete set and maximizing a penalty function; 4, performing single-step correction on the initial estimation values of the current frequency and complex amplitude by adopting a Newton method; step 5, circularly correcting the frequency and the complex amplitude value of the signal; step 6, updating the complex amplitude value through a least square method; 7, detecting whether the residual observation energy is smaller than a threshold value or not; and step 8, outputting a target estimation result set. The method has the advantages of being high in resolution, high in estimation accuracy, capable of reducing the influence of lattice point mismatch and low in complexity.

Description

Technical Field

[0001] The invention belongs to the field of array signal processing, and in particular relates to a two-dimensional super-resolution detection method for solving an arrival direction problem using an inverted echelon array. Background Art

[0002] Locating multiple targets, or determining the azimuth and elevation angles of each signal source, is an important problem in the field of signal processing.

[0003] The early traditional beamforming algorithm is a powerful low-complexity algorithm for estimating the target position. It determines the target position by performing a Fourier transform on the received array data. This method is a simple spatial expansion process in the traditional Fourier spectrum estimation algorithm in the time domain, that is, the algorithm uses the data received by each array element in the spatial domain to replace the time domain data in the traditional time domain processing, but its resolution is relatively low.

[0004] Based on the beamforming algorithm, scholars have proposed some high-resolution algorithms based on subspace, such as pencil-MUSIC2 and 2D unitary ESPRIT. The principle of the MUSIC algorithm is as follows: (1) obtain the data covariance matrix based on the received data; (2) perform eigendecomposition on the data covariance matrix; (3) determine the number of signal sources based on the eigenvalues ​​of the data covariance matrix; (4) determine the signal subspace and noise subspace; (5) perform spectrum peak search based on the signal parameter range; (6) find the angle corresponding to the maximum point, which is the signal incident direction.

[0005] A common feature of this type of algorithm is that the received data is divided into two mutually orthogonal subspaces through mathematical decomposition of the array received data (such as eigendecomposition, singular value decomposition and QR decomposition, etc.): one is the signal subspace consistent with the array flow space of the signal source, and the other is the noise subspace orthogonal to the signal subspace. Subspace decomposition algorithms use the orthogonal characteristics of the two subspaces to construct a "pin-volt" spatial spectrum peak, thereby greatly improving the resolution of the algorithm. However, the performance of these subspace-based high-resolution algorithms is very sensitive to the number of signal sources, and the algorithm performance will deteriorate rapidly under low signal-to-noise ratio conditions.

[0006] Later, scholars took advantage of the sparsity of signal sources in the spatial domain and proposed some methods based on compressed sensing. Due to its strong theoretical basis and high resolution, the method based on compressed sensing has become an attractive algorithm. The traditional method based on compressed sensing builds a discrete network and assumes that all estimated positions fall precisely on the grid points. However, this assumption is unrealistic, and when the estimated position does not fall on the grid points, a grid mismatch problem occurs. In order to overcome this grid matching problem, it is necessary to do a good job of grid division, but this will greatly increase the computational complexity.

[0007] Another signal estimation method is to regard the target's azimuth as a parameter that can be continuously valued, which is a method based on atomic norm minimization. However, the algorithm based on the atomic norm involves solving a positive semidefinite programming problem, so the computational complexity is very high.

[0008] For the problem of estimating the two-dimensional arrival direction based on a uniform circular array, the early traditional beamforming algorithm requires precise design of the position and layout of multiple sensors. Improper layout may affect the beamforming effect. In addition, only one Fourier transform operation is performed on the array data, so the resolution is relatively low. The performance of the high-resolution algorithm based on subspace is very sensitive to the number of signal sources, and the algorithm performance will deteriorate rapidly under low signal-to-noise ratio conditions. The method based on compressed sensing will have grid mismatch problems because it builds a discrete network. The method based on atomic norm minimization has high computational complexity because it involves matrix inversion. Summary of the invention

[0009] The purpose of the present invention is to provide a two-dimensional super-resolution detection method for solving the two-dimensional arrival direction estimation problem based on an inverted echelon array, so as to better estimate the target azimuth, elevation angle and complex amplitude vector, and has the advantages of high resolution, high estimation accuracy, ability to reduce the influence of grid mismatch and low complexity.

[0010] The technical solution of the present invention is:

[0011] A two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array comprises the following steps:

[0012] Step 1, initialization: set the stop threshold τ, the maximum number of iterations I, and the number of single-step Newton corrections R s and the number of cyclic Newton corrections R c The initial parameters of the algorithm;

[0013] Step 2: Input array data: Input the array data Y received by the equally spaced circular array, straighten it to y, i = 0, P 0 ={}, when And when i≤I, repeat the following steps: i=i+1;

[0014] Step 3, detection phase: use a finite discrete set instead of an infinite set to obtain an initial estimate of the line spectrum frequency w by maximizing the penalty function;

[0015] Step 4: Single-step correction: Use Newton's method (using the first-order derivative and second-order derivative of the penalty function) to perform single-step correction on the initial estimate of the current frequency and complex amplitude in step 3, and repeat it for R times. s , get the updated target information;

[0016] Step 5, cyclic correction: re-correct the frequency and complex amplitude value of the signal in step 4, fix the information of other signals, perform a Newton correction on the frequency and complex amplitude value of each signal in turn, and obtain the updated target information. The number of cycles is R c ;

[0017] Step 6, least squares correction: Update the complex amplitude value through the least squares method to ensure that the residual energy is the minimum possible value of the current set of estimated frequencies, and let P m Represents a new set of parameters;

[0018] Step 7, residual energy threshold judgment: remove the effective target estimation result signal data from the initial signal, and detect whether the residual observation energy is ≤ the threshold τ;

[0019] Step 8: Program stop: Stop calculation after reaching the set conditions and output the target estimation result set.

[0020] In step 3, w is restricted to a discrete set (Ω θ ,Ω φ ), (γ θ ,γ φ ) is the oversampling rate, and w is defined c =[θ c ,φ c ] H = argmax G Y (w) is the output of this stage, and the corresponding estimate of the complex amplitude x(t) is At this time, after detection

[0021] In step 4, let Represents the current estimate, then the framework of the frequency-corrected Newton orthogonal matching pursuit algorithm is:

[0022]

[0023] in:

[0024]

[0025] and

[0026]

[0027] yes right In terms of the gradient and Hansen matrix, Indicates partial derivative. At this time, after correction

[0028] Further, in step 4,

[0029]

[0030] represents the derivative operation;

[0031] When the penalty function is a locally concave function, the above update rule is used to maximize

[0032]

[0033] It is the straightening of matrix A, each column of A is

[0034]

[0035] ☉ represents the Hadamard product, which means multiplying corresponding elements of two matrices of equal dimension.

[0036] In step 4, the method for estimating multiple line spectrum frequencies is to perform Newton coordinate descent on the residual energy using a single signal correction procedure. For the estimated L signals, let P = {(x l ,w l ),l=1,…,L} represents the estimated value of the detected target, then the straightened matrix of the residual observation corresponding to this estimated value is defined as y r (P), then:

[0037]

[0038] The process of correcting the lth signal is as follows: r (P\{x l ,ω l}) is considered as the initial observation y, and then the single-step Newton correction step is used to update (x l ,ω l ), the next step based on coordinate descent involves adjusting all detected frequencies w l .

[0039] In step 4, the correction acceptance criteria are: the correction step is accepted if the correction will result in Strict improvement of

[0040] In step 5, in each Newton correction, Treat it as observation Y, and then apply the single-step correction step. In step 6, the least squares correction formula is:

[0041] x l (t)(t=1,...,T,l=1,...,m)

[0042] ·) + It means taking the pseudo-inverse and projecting the received signal onto the subspace spanned by the estimated frequency, so as to obtain the updated frequency and amplitude values.

[0043] In step 7, estimated signal data is generated based on the frequency and complex amplitude values ​​in all current target information, and these estimated signal data are removed from the initial data of the last observation, and the remaining residual observation energy is used as the new initial observation data y.

[0044] In step 8, the conditions for stopping repetition are: a. the number of repetitions has reached the maximum number of iterations preset by the algorithm; b. the algorithm determines that the residual observation energy does not exceed the corresponding threshold.

[0045] The beneficial effects of the present invention are:

[0046] 1. The present invention uses the Newton gradient descent method to perform single-step and cyclic corrections on the grid point estimation results, which effectively alleviates the grid point mismatch problem and improves the frequency estimation accuracy;

[0047] 2. The present invention utilizes the advantages of fast Fourier transform and loop iteration to completely avoid the operation of matrix inversion, greatly reducing the computational complexity and facilitating chip implementation;

[0048] 3. In order to solve the problem of mutual interference between targets in space, the present invention adopts a cyclic correction step to correct all target results in turn, which effectively reduces the mutual influence and further improves the estimation accuracy;

[0049] 4. The present invention can use multiple snapshots of observation data to improve estimation accuracy. BRIEF DESCRIPTION OF THE DRAWINGS

[0050] Figure 1 It is a model diagram of the two-dimensional arrival direction estimation problem based on the inverted echelon array in the present invention;

[0051] Figure 2 It is a workflow diagram of a two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted ladder array of the present invention. DETAILED DESCRIPTION

[0052] The present invention will be further described below in conjunction with specific embodiments and accompanying drawings:

[0053] like Figure 1 As shown in the figure, a problem of estimating the arrival direction of an inverted echelon array is presented. It is assumed that there is an inverted echelon array with M layers and N array elements arranged at equal intervals in each layer. The radius of each layer is r m , which receives the observation data of K independent narrowband signals in the far field. Let λ represent the signal wavelength, (θ k ,φ k ) represents the pitch angle and azimuth angle of the kth target, θ k ∈[0,π / 2],φ k ∈[0,2π]. The first array element is placed on the positive half axis of the bottom layer x-axis, and the angle between the nth array element and the positive half axis of the x-axis is β n , the surface position vector is p n =(r m cosβ n ,r m sinβ n ,0), the unit vector of the arrival direction is q k =(sinθ k cosφ k ,sinθ k sinφ k ,cosθ k ).

[0054] The interval between the angles formed by any two array elements and the center O is 2π / M, so The delay difference between the center O and the nth array element is:

[0055]

[0056] Here, · represents the dot product operation. To perform a dot product operation on two vectors, we multiply the elements at corresponding positions of the two vectors one by one and then sum them. (·) T Represents a transpose operation;

[0057] Combining c = fλ, c is the speed of sound, f is the frequency, the phase difference between the center O and the nth array element is:

[0058]

[0059] Since most signals are transmitted from far away, we can make the following assumption: for a target, the centers of different layers are consistent, that is, for all layers, (θ k ,φ k ,β n ) is the same. It is worth noting that due to limited computing resources, in practice, not all array elements in each layer are used for signal estimation. In addition, some array elements may not be able to receive signals due to occlusion by certain structures.

[0060] Here we assume that each layer has only N B Sensors are used for beamforming. For example, if the signal is incident on the nth sensor, then the left and right N sensors of each layer are B / 2 sensors are used for estimation.

[0061] At the tth beat (taking the center as the reference point), the observation can be modeled as:

[0062]

[0063] in Expressed by all M layers, each layer N B The observation data collected by the sensors, x k (t) represents the complex magnitude vector, which has the same magnitude in each layer. The mth column

[0064]

[0065] is the array flow pattern vector; is additive white Gaussian noise, each element is an independent and identically distributed complex Gaussian random variable with mean 0 and variance σ 2 ; x(t) = [x 1 (t),x 2 (t),…,x K (t)] T .

[0066] for Figure 1 The problem model shown, assuming there is only one goal, the model is simplified to:

[0067]

[0068] in The mth column

[0069]

[0070] For simplicity, normalize it:

[0071]

[0072] Straighten the left and right sides of the model at the same time, then:

[0073] y(t)=a(w)x(t)+n(t

[0074] in y=vec(Y),a(w)=vec(A(w)),n(t)=vec(N(t))

[0075] The maximum likelihood estimation of x and w can be obtained by minimizing the residual energy To calculate, here represents the square of the F norm, which means summing the squares of the moduli of all elements of the matrix. This process is equivalent to maximizing S(x,w).

[0076]

[0077] Here Re represents the real part operation, || represents the modulus, ||·|| represents the two-norm operation, (·) H represents conjugate transpose, and ∑(·) represents summation.

[0078] It is difficult to optimize this function directly, therefore, the present invention takes two steps:

[0079] In the first step, the detection phase, a rough estimate is determined by oversampling spectrum peak search;

[0080] The second step, the correction phase, iteratively updates the amplitude and frequency of the estimated signal.

[0081] For any given w, the gain of maximizing S(x(t),w) is Will Substituting into S(x,ω), we can obtain the generalized likelihood ratio detection estimate for the frequency, which is equivalent to the solution to the following optimization problem:

[0082]

[0083] Here arg is the abbreviation of argument, and argmax is the value of the variable when the following formula reaches the maximum value. y(t) is the GLRT penalty function based on the t-th snapshot.

[0084] Based on the above equivalent estimation of the problem model, in order to find the global optimal solution to the line spectrum estimation problem based on the inverted terrace array, the present invention provides a Newton orthogonal matching pursuit algorithm for detection and correction, such as Figure 2 As shown, the specific steps include:

[0085] Step 1, initialization: set the stop threshold τ, the maximum number of iterations I, and the single step R s and the number of cyclic Newton corrections R c The initial parameters of the algorithm;

[0086] Step 2: Input array data: Input the array data Y received by the equally spaced circular array, straighten it to y, i = 0, P 0 ={}, when And when i≤I, repeat the following steps: i=i+1;

[0087] Step 3: Detection phase: By restricting w to a discrete set (Ω θ ,Ω φ ), the present invention can obtain a rough estimate of the line spectrum frequency w, where (γ θ ,γ φ ) is the oversampling rate. At the same time, define w c =[θ c ,φ c ] H =argmaxG Y 9w is the output of this stage, and the corresponding estimate of the complex amplitude x(t) is At this time, after detection

[0088] Step 4, single-step correction stage: let represents the current estimated value, then the Newton process of frequency correction is,

[0089]

[0090] in

[0091]

[0092] and

[0093]

[0094] yes right In terms of the gradient and Hansen matrix, It means to find partial derivative;

[0095] At this time, after correction

[0096] Specifically,

[0097]

[0098]

[0099] Represents a derivative operation.

[0100] When the penalty function is a local concave function, the present invention adopts the above update rule to maximize

[0101] It is the straightening of matrix A, and each column of A is

[0102]

[0103] ☉ represents the Hadamard product, which means multiplying corresponding elements of two matrices of equal dimension.

[0104] The above gives the framework of the Newton orthogonal matching pursuit algorithm assuming that there is only one target. Next, for the L signals that have been estimated, let P = {9x l ,w l ), l = 1, ..., L} represents the estimated value of the detected target, then the straightened matrix of the residual observation corresponding to this estimated value is defined as y r (P)

[0105]

[0106] The method for estimating multiple line frequencies is to perform Newton coordinate descent on the residual energy using a single signal correction procedure. One step of this coordinate descent involves adjusting all detected frequencies w l .

[0107] The process of correcting the lth signal is as follows: r (P\{x l ,ω l}) is considered as the initial observation y, and then the single-step Newton correction step is used to update (x l ,ω l ). Correction acceptance criteria: Correction steps are accepted if the correction leads to Strict improvement, that is, By doing so, the present invention can ensure that the correction steps accepted will definitely result in a reduction of the residual energy.

[0108] Step 5: Cyclic correction phase: The optimization step can be implemented multiple times. Single-step correction refers to using the gradient of the penalty function and the Hansen matrix to correct the estimation result previously selected roughly in the dictionary matrix, the number of times is R s Cyclic correction can be regarded as a single-step correction of each signal in all detected signals, the number of times is Rc . Let R s and R c Respectively represent the number of single-step correction and cyclic correction.

[0109] For all P" m Specifically, the information of other signals is fixed, and the frequency and complex amplitude value of each signal are subjected to a Newton correction in turn. is regarded as observation Y, and then a single-step correction step is applied to obtain the updated target information.

[0110] Step 6, least squares correction: Update the complex amplitude value by least squares method:

[0111] x l (t)(t=1,...,T,l=1,...,m)

[0112] ·) + It means taking the pseudo-inverse and projecting the received signal onto the subspace spanned by the estimated frequency, so that the updated frequency and amplitude values ​​are obtained, and at the same time, let P m The present invention updates the amplitude of the signal to ensure that the residual energy is the minimum possible value for the current set of estimated frequencies.

[0113] Step 7, residual energy threshold judgment: Eliminate the effective target estimation result signal data from the initial signal, and detect whether the residual observation energy is ≤ the threshold τ. Specifically, based on the frequency and complex amplitude values ​​in all the current target information, generate the estimated signal data, eliminate these estimated signal data from the initial data of the last observation, and use the remaining residual observation energy as the new initial observation data y.

[0114] Step 8, program stop: stop the calculation after reaching the set conditions, and output the target estimation result set. Specifically, the conditions for repeated stopping are: a. the number of repetitions has reached the maximum number of iterations preset by the algorithm; b. the algorithm determines that the residual observation energy does not exceed the corresponding threshold.

[0115] In summary, the present invention first detects an estimate of a frequency in a discrete set pair by maximizing a penalty function, and then updates its corresponding amplitude estimate. Then, the estimate is improved using information from the first and second order derivatives of the penalty function. Next, the estimate of each previously detected target is further improved using information from all other previously detected targets, one at a time, which is critical to the convergence and accuracy of the algorithm. Finally, the amplitude is updated by the least squares method.

[0116] The advantages of the present invention are: a two-dimensional super-resolution detection method is proposed to solve the problem of arrival direction estimation based on an inverted ladder array. The method utilizes the advantages of fast Fourier transform and iterative loop, avoids matrix inversion, has low computational complexity, and is easy to implement in hardware. At the same time, the method uses Newton single-step and cyclic correction feedback strategies to reduce the impact of grid mismatch problems, reduce mutual influence between targets, have high estimation accuracy, and can use multi-shot data to achieve near-optimal estimation performance.

[0117] It should be noted that the above embodiments are merely preferred embodiments of the present invention and are not intended to limit the protection scope of the present invention. All equivalent modifications made on the basis of the above embodiments belong to the protection scope of the present invention.

Claims

1. A two-dimensional super-resolution detection method using an inverted echelon array to solve the arrival direction problem, characterized in that: The steps include: Step 1, initialization: set the stop threshold τ, the maximum number of iterations I, and the number of single-step Newton corrections R s and the number of cyclic Newton corrections R c The initial parameters of the algorithm; Step 2, array data input: input the array data Y received by the equally spaced circular array, straighten it to y, i = 0, P0 = {}, when And when i≤I, repeat the following steps: i=i+1; Step 3, detection phase: use a finite discrete set instead of an infinite set to obtain an initial estimate of the line spectrum frequency w by maximizing the penalty function; Step 4: Single-step correction: Use Newton's method to perform single-step correction on the initial estimated values ​​of the current frequency and complex amplitude in step 3, and repeat the correction number R. s , get the updated target information; Step 5, cyclic correction: re-correct the frequency and complex amplitude value of the signal in step 4, fix the information of other signals, perform a Newton correction on the frequency and complex amplitude value of each signal in turn, and obtain the updated target information. The number of cycles is R c ; Step 6, least squares correction: Update the complex amplitude value through the least squares method to ensure that the residual energy is the minimum possible value of the current set of estimated frequencies, and let P m Represents a new set of parameters; Step 7, residual energy threshold judgment: remove the effective target estimation result signal data from the initial signal, and detect whether the residual observation energy is ≤ the threshold τ; Step 8: Program stop: Stop calculation after reaching the set conditions and output the target estimation result set.

2. The two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array according to claim 1, characterized in that: In step 3, w is restricted to a discrete set (Ω θ ,Ω φ ), (γ θ ,γ φ ) is the oversampling rate, and w is defined c =[θ c ,φ c ] H = argmax G Y (w) is the output of this stage, and the corresponding estimate of the complex amplitude x(t) is At this time, after detection 3. The two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array according to claim 2, characterized in that: In step 4, let Represents the current estimate, then the framework of the frequency-corrected Newton orthogonal matching pursuit algorithm is: in: and yes right In terms of the gradient and Hansen matrix, Indicates partial derivative. At this time, after correction 4. The two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array according to claim 3 is characterized by: In step 4, represents the derivative operation; When the penalty function is a locally concave function, the above update rule is used to maximize It is the straightening of matrix A, and each column of A is ☉ represents the Hadamard product, which means multiplying corresponding elements of two matrices of equal dimension.

5. The two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array according to claim 4, characterized in that: In step 4, the method for estimating multiple line spectrum frequencies is to perform Newton coordinate descent on the residual energy using a single signal correction procedure. For the estimated L signals, let P = {(x l ,w l ),l=1,…,L} represents the estimated value of the detected target, then the straightened matrix of the residual observation corresponding to this estimated value is defined as y r (P), then: The process of correcting the lth signal is as follows: r (P\{x l ,ω l }) is considered as the initial observation y, and then the single-step Newton correction step is used to update (x l ,ω l ), the next step based on coordinate descent involves adjusting all detected frequencies w l .

6. The two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array according to claim 5, characterized in that: In step 4, the correction acceptance criteria are: the correction step is accepted if the correction will result in Strict improvement of 7. The two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array according to claim 6, characterized in that: In step 5, in each Newton correction, is considered as observation Y, and then a single correction step is applied.

8. The two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array according to claim 7, characterized in that: In step 6, the least squares correction formula is: x(t)(t=1,……,T,l=1,……,m) ·) + It means taking the pseudo-inverse and projecting the received signal onto the subspace spanned by the estimated frequency, so as to obtain the updated frequency and amplitude values.

9. The two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array according to claim 8, characterized in that: In step 7, estimated signal data is generated based on the frequency and complex amplitude values ​​in all current target information, and these estimated signal data are removed from the initial data of the last observation, and the remaining residual observation energy is used as the new initial observation data y.

10. A two-dimensional super-resolution detection method for solving the arrival direction problem using an inverted echelon array according to any one of claims 1 to 9, characterized in that: In step 8, the conditions for stopping repetition are: a. the number of repetitions has reached the maximum number of iterations preset by the algorithm; b. the algorithm determines that the residual observation energy does not exceed the corresponding threshold.

Images