A Direction of Arrival Estimation Method Based on Array Signal Processing

Abstract

The present invention discloses a method for estimating a direction of arrival based on array signal processing. The method proposes a cyclostationary signal orthogonal matching pursuit algorithm based on rotationally invariant subspace information for indoor positioning scenarios. A cyclostationary algorithm is introduced as the initial correlation function to combine the cyclostationarity on the time scale with the array signal on the space scale, thereby effectively solving the problem that the orthogonal matching pursuit algorithm is sensitive to noise and has poor noise resistance; rotationally invariant subspace information is used to break through the Rayleigh limit of the orthogonal matching pursuit algorithm, effectively improving the multi-object positioning resolution and reducing the computational complexity; an orthogonal matching pursuit algorithm based on compressed sensing theory is used to further improve the algorithm's anti-interference ability through sparse signal reconstruction. The requirements of high precision and low latency for indoor positioning scenarios are met.

Description

Technical Field

[0001] The invention belongs to the technical field of direction of arrival positioning, and in particular relates to a direction of arrival estimation method based on array signal processing. Background Art

[0002] As an important research direction in the field of signal processing, array signal processing has received widespread attention and rapid development in recent years. This technology obtains discrete observation data of signals by placing a signal receiving sensor array in a certain flow pattern in three-dimensional space, and uses the differences in data received by different sensors to design corresponding algorithms to obtain various parameters of the signal. In recent years, array signal processing has been widely used in many civil or military fields such as sonar communication, radio astronomy, biomedicine, positioning and navigation.

[0003] The research on array signal processing mainly includes two aspects: adaptive spatial filtering and spatial spectrum estimation. Among them, spatial spectrum estimation is also called Direction Of Arrival (DOA) estimation. The DOA estimation problem is to determine the spatial positions of multiple signals simultaneously located in a certain area of ​​space based on the information of each signal received by the array reference antenna, that is, the azimuth angle of each signal relative to the array reference antenna. The resolution of DOA estimation is related to the length of the antenna array. When the length of the antenna array is determined, its corresponding resolution is also determined, which is called the Rayleigh limit. For the Rayleigh limit, there are methods that can break through the Rayleigh limit. Such methods are called super-Rayleigh limit methods or super-resolution methods.

[0004] Early DOA estimation research only studied the one-dimensional arrival angle of the spatial incoming signal, and the receiving array usually used a uniform linear array. The first proposed method was the conventional beamforming method. Subsequently, many estimation methods that can break through the "Rayleigh limit" and achieve higher resolution have been proposed by many people. Among them, the multiple signal classification (MUSIC) algorithm is the most famous. In 1989, the rotation invariant technology (ESPRIT) algorithm was proposed by Roy R et al. In the 1980s, subspace fitting algorithms such as the maximum likelihood (ML) algorithm and the subspace fitting (SF) algorithm appeared. Compared with algorithms such as the MUSIC algorithm and the ESPRIT algorithm, these algorithms estimate the angle values ​​more accurately, but the computational complexity is also greater.

[0005] In recent years, the development of compressed sensing (CS) theory has also provided new ideas for DOA estimation in array signal processing. CS theory breaks through the limitations of the Nyquist theorem. When the signal is sampled at a sampling frequency lower than the Nyquist frequency, it can completely restore the original signal from the sampled signal with a small number of dimensions. When CS theory is applied to spatial spectrum estimation, the received signal can be regarded as a signal sampled at a frequency lower than the Nyquist frequency in the entire spatial domain, and then the arrival angle of the source can be found by solving the optimal solution. Compared with traditional DOA estimation methods, DOA estimation based on compressed sensing theory (CS-DOA) has many advantages, including low sampling rate, shorter sampling time, strong anti-interference ability, etc. At present, CS-DOA estimation is still in the stage of continuous exploration, but its theoretical research value and practical engineering significance are immeasurable. Summary of the invention

[0006] The present invention aims to provide a method for estimating the direction of arrival based on array signal processing, wherein incoming wave sampling data is obtained through an antenna array, a cyclic correlation matrix is ​​calculated by a cyclostationary algorithm, and then a corresponding signal subspace is obtained by using singular value decomposition, a sparse reconstruction model is constructed, and then a dictionary vector corresponding to each iterative residual is obtained based on the rotationally invariant subspace information, the dictionary vector is matched and added to a support set to calculate and update the residual, and it is continuously iterated until all incoming wave directions are calculated, so as to solve the technical problem that the requirements for accuracy and delay of indoor positioning cannot be met.

[0007] In order to solve the above technical problems, the specific technical solutions of the present invention are as follows:

[0008] A method for estimating a direction of arrival based on array signal processing comprises the following steps:

[0009] Step 1: Obtain incoming wave sampling data through the antenna array;

[0010] Step 2, calculating the cyclic cross-correlation matrix by a cyclostationary algorithm;

[0011] Step 3: Use singular value decomposition to obtain the corresponding signal subspace and construct a sparse reconstruction model;

[0012] Step 4: Obtain the dictionary vector corresponding to each iteration residual based on the rotation invariant subspace information;

[0013] Step 5: Match the dictionary vectors and add them to the support set to calculate the updated residual;

[0014] Step 6: Continue iterating until all incoming wave directions are calculated. Repeat steps 4 and 5, and continue iterating until all incoming wave directions are calculated.

[0015] Further, step 1 specifically includes the following steps: using a uniform linear array to receive the signal, for a uniform linear array with N array elements spaced d apart, f 0 is the signal frequency, λ 0 is the signal wavelength; when there are M signal sources with a cyclic frequency of α incident, where M≤N, the signal s received by the kth array element at the tth snapshot is obtained k (t), the signal wave direction θ i Satisfies -π / 2≤θ i ≤π / 2, where i=1,2,…,M, indicating that there are M signal sources in total; the array sampled data S(t)=A(θ)X(t)+N(t) is obtained, t=1,2,…,U, where A(θ) represents the received signal angle matrix, X(t) represents the basic signal matrix, N(t) represents the noise matrix, t represents the tth snapshot data of the array, and it is assumed that a total of U snapshot data are received.

[0016] Furthermore, the calculation method of step 2 includes the following steps:

[0017] Step 2.1: Calculate the time delay τ between the kth array element and the i-th signal source ki ;

[0018] Step 2.2: Based on the definition of cyclic cross-correlation and the array element receiving data, the array element s can be obtained. p and q Circular cross-correlation function of the received data

[0019] Step 2.3: Where p,q=1,2,.,N, solve the cyclic cross-correlation function for N different array elements to obtain the cyclic cross-correlation matrix of the array received signal Where p and q represent the pth and qth array elements.

[0020] Furthermore, in step 2.1,

[0021] In step 2.2,

[0022] in represents the cyclic autocorrelation of the basic signal of the ith signal source, τ pi represents the delay between the pth array element and the i-th signal source, τ qi represents the time delay between the qth array element and the i-th signal source, and for an equidistant uniform linear array, the following formula holds true:

[0023]

[0024] Since the noise cycle frequency is not α, its correlation coefficient at the α frequency is 0, and the cyclic cross-correlation function is:

[0025]

[0026] where s p (t) represents the data sampled from the pth array, s q (t) represents the data of the qth array sample;

[0027] The cyclic cross-correlation matrix of the array received signal in step 2.3 is calculated as:

[0028]

[0029] Further, step 3 includes the following steps:

[0030] For the circulant cross-correlation matrix obtained in step 2 Perform singular value decomposition to obtain unitary matrices U and V; arrange the singular values ​​in order from large to small, and divide the unitary matrices U and V into two parts: signal subspace and noise subspace according to the size of the singular values; span the signal subspace U from the eigenvectors of the M largest singular values S and V S ; Finally, the new signal output equation R is obtained S .

[0031] Furthermore, in step 3, the cyclic cross-correlation matrix The i-th column of is reconstructed like the sparse recovery model: Where i = 1, 2, ..., N; and r i Both said The i-th column of is the dictionary matrix, b i It is a dictionary Next i The sparse coefficient vector of ; represents the noise value, e i is an N×1 column vector, the i-th element is 1 and the rest are 0;

[0032] according to Change it to matrix form as follows: Where B = [b 1 ,b 2 ,...,b N ] and b i (i=1,2,..,N) have the same sparse structure, that is, each non-zero element of the ideal value should appear in the same row of B, and the number of rows filled with 0 is the same, I N is the N×N identity matrix;

[0033] in is an N×N matrix; The singular value decomposition of is expressed as: Where U and V are both unitary matrices, and Σ consists of singular values; therefore There are M large singular values; arrange the singular values ​​in order from large to small: λ 1 ≥λ 2 ≥...≥λ M ≥λ M+1 ≥...≥λ N ; The last NM singular values ​​are caused by noise and interference; according to the size of the singular values, the unitary matrices U and V are divided into two parts: where Σ S It is composed of the first M largest singular values, Σ N It is composed of the remaining singular values;

[0034] Among them U S and V S represents the signal subspace spanned by the eigenvectors of the M largest singular values, and the Σ diagonal The remaining NM singular values ​​of constitute the noise subspace U N and V S ; A new signal output equation can be obtained:

[0035]

[0036] in is the noise interference error term; ΔR is given by and The error caused by the approximation of is a component of the signal subspace.

[0037] Furthermore, the calculation method of step 4 includes the following steps:

[0038] Step 4.1: Initialize the state support set S 0 =φ, residual C 0 =R S , the sparse representation matrix

[0039] Step 4.2: Let this be the kth iteration: the support set is S k , the residual is C k , the sparse representation matrix is The first P most significant solutions of the kth iteration are θ i ,i=1,2,...,P; get A(θ)·T=C k , where T is the rotation matrix; let

[0040]

[0041] Where A(θ) represents the received signal angle matrix, A 1 (θ) represents the first N-1 rows of A(θ), A 2 (θ) represents the last N-1 rows of A(θ), Represents C k The first N-1 rows of the matrix, Respectively represent C k The last N-1 rows of the matrix;

[0042] Let the transition matrix Ψ, where A 2 (θ) = A 1 (θ)·Ψ; calculated Due to T -1 ΨT is similar to Ψ and Ψ is a diagonal matrix, so Ψ can be solved; thus the dictionary vector α of the kth iteration is calculated i ,i=1,2,...,P.

[0043] Further, step 5 includes the following steps: let the index of the vector selected in the kth iteration be i k ,Pick will i k Add to S k , update the support set and sparse representation matrix; thereby updating the residual iteration equation of orthogonal matching pursuit and obtaining a new residual solution C k+1 .

[0044] Furthermore, in step 5, the support set and the sparse representation matrix are updated: S k+1 =S k ∪{i k},

[0045] The residual iteration equation of orthogonal matching pursuit is expressed as:

[0046]

[0047]

[0048] The subscript k represents the kth iteration, B represents the unknown number of the residual iteration equation, and B k represents the least squares solution of the kth iteration, C k represents the residual signal matrix of the kth iteration, R S represents the signal output equation obtained in step 3, represents the new dictionary matrix consisting of the k atoms selected at the kth iteration, express The pseudo-inverse matrix, F refers to the two norm; from the above formula, the corresponding residual Ck+1 , and iterate.

[0049] The method for estimating direction of arrival based on array signal processing of the present invention has the following advantages:

[0050] The present invention is based on a cyclostationary algorithm, rotationally invariant subspace information and an orthogonal matching pursuit algorithm to obtain a wave direction estimation method based on array signal processing. The cyclostationary algorithm is used to combine the cyclostationarity on the time scale with the array signal on the space scale, so that the anti-noise performance of the algorithm is enhanced. The method breaks through the Rayleigh limit of the orthogonal matching pursuit algorithm by combining the rotationally invariant subspace information with the orthogonal matching pursuit algorithm based on the compressed sensing theory, effectively improves the multi-object positioning resolution and reduces the computational complexity, and further improves the anti-interference performance of the algorithm. The requirements of high precision and low latency for indoor positioning scenarios are met, and the user experience is improved. BRIEF DESCRIPTION OF THE DRAWINGS

[0051] Figure 1 A schematic diagram of a scenario of a method for estimating a direction of arrival based on array signal processing according to the present invention;

[0052] Figure 2 It is a flow chart of the method for estimating the direction of arrival based on array signal processing of the present invention;

[0053] Figure 3 The flowchart of the Cyclic-orthogonal matching pursuit algorithm based on the direction of arrival of array signal processing of the present invention;

[0054] Figure 4 It is a flow chart of the simplified dictionary orthogonal matching pursuit algorithm based on rotation-invariant subspace information of the direction of arrival based on array signal processing of the present invention. DETAILED DESCRIPTION

[0055] In order to better understand the purpose, structure and function of the present invention, a method for estimating a direction of arrival based on array signal processing of the present invention is further described in detail below with reference to the accompanying drawings.

[0056] like Figure 1 As shown, the research of the present invention is based on indoor positioning scenarios. The transmitter bound to the object to be positioned is used to send signals, multiple antenna arrays at different positions are used to sample received signals, and the processor is used to operate according to the sampled signals to execute the method described, and after the obtained angle is input into the PC processing end, high-precision positioning information is obtained according to multiple angle data.

[0057] like Figure 2As shown, the present invention designs a method for Bluetooth wave arrival direction based on array signal processing, through which the signal is processed to obtain the corresponding wave arrival direction, which is used to achieve high-precision and low-latency indoor positioning. In practical application, it specifically includes the following steps.

[0058] Step 1: Obtain incoming wave sampling data through the antenna array. We use a uniform linear array to receive the signal. For a uniform linear array with N array elements separated by d, when only one signal source is incident and noise is ignored, the signal received by the first array element is: Where A(t) represents the amplitude of the received signal, e is the base of the natural logarithm function, j is an imaginary number, and t is the corresponding time. The delay between adjacent array elements is Where c is the speed of light and θ is the direction of the corresponding wave. In a narrowband signal, the second array element receives the signal:

[0059]

[0060] Among them, f 0 is the signal frequency, λ 0 is the signal wavelength.

[0061] When there are M signal sources with a cyclic frequency of α incident (M≤N), the signal received by the kth array element at the tth snapshot is:

[0062] in A represents the signal received by the kth array element in the direction of the i-th incoming wave. i (t) represents the amplitude of the received signal, represents the time delay of the kth array element when receiving the ith incoming wave direction, n k (t) is the noise received by the kth array element. Its incoming wave direction θ i Satisfies -π / 2≤θ i ≤π / 2, where i=1,2,…,M. Let Then the data sampled by the array at the tth snapshot can be expressed as:

[0063]

[0064] That is: S(t)=A(θ)X(t)+N(t), t=1,2,......,U. Where U is the number of snapshots, A(θ)=[a(θ 1 ),...,a(θ M )],a(θ i ) is the direction vector of the i-th source.

[0065] In step 2 and step 3, the Bluetooth direction of arrival Cyclic-orthogonal matching pursuit algorithm based on array signal processing designed by the present invention is described, such as Figure 3 shown.

[0066] Step 2: Calculate the cyclic correlation matrix using the cyclostationary algorithm. This algorithm uses the Cyclic algorithm for calculation. Assume that the time delay τ ki , where k represents the kth array element, i represents the ith signal source, According to the definition of cyclic cross-correlation and the array element receiving data, the array element s can be obtained p and q The cyclic cross-correlation function of the received data is

[0067]

[0068] where p = 1, 2, ..., N, q = 1, 2, ..., N, represents the cyclic autocorrelation of the basic signal of the ith signal source, τ pi represents the delay between the pth array element and the i-th signal source, τ qi represents the time delay between the qth array element and the i-th signal source. And for equidistant uniform linear arrays, the following formula holds true:

[0069]

[0070] Since the noise cycle frequency is not α, its correlation coefficient at the α frequency is 0, and the cyclic cross-correlation function is:

[0071]

[0072] where s p (t) and s q (t) represent the data of the pth and qth array samples respectively.

[0073] because Where p,q=1,2,.,N, then the cyclic correlation matrix of the array receiving signal can be calculated as:

[0074]

[0075] Step 3: Use singular value decomposition to obtain the corresponding signal subspace, build a sparse reconstruction model, and circulate the correlation matrix The i-th column of is reconstructed like the sparse recovery model: Where i = 1, 2, ..., N; and r i Both said The i-th column of is the dictionary matrix, b i It is a dictionary Nexti The sparse coefficient vector of ; represents the noise value, e i is an N×1 column vector, the i-th element is 1 and the rest are 0;

[0076] according to Change it to matrix form as follows: Where B = [b 1 ,b 2 ,...,b N ] and b i (i=1,2,..,N) have the same sparse structure, that is, each non-zero element of the ideal value should appear in the same row of B, and the number of rows filled with 0 is the same, I N is the N×N identity matrix;

[0077] in is an N×N matrix; The singular value decomposition of can be expressed as: Where U and V are both unitary matrices, and Σ consists of singular values; therefore There are M large singular values; arrange the singular values ​​in order from large to small: λ 1 ≥λ 2 ≥...≥λ M ≥λ M+1 ≥...≥λ N ; The last NM singular values ​​are caused by noise and interference; according to the size of the singular values, the unitary matrices U and V are divided into two parts: where Σ S ,Σ N They are respectively composed of the first M largest singular values ​​and the remaining singular values.

[0078] Among them U S (V S ) represents the signal subspace spanned by the eigenvectors of the M largest singular values, and the Σ diagonal The remaining NM singular values ​​of constitute the noise subspace U N (V S ); Combining the relationship shown in the above two equations, a new signal output equation can be obtained:

[0079]

[0080] in is the noise interference error term; ΔR is given by and The error caused by the approximation of is a component of the signal subspace.

[0081] In step 4, step 5 and step 6, the simplified dictionary orthogonal matching pursuit algorithm based on rotation-invariant subspace information of Bluetooth direction of arrival based on array signal processing designed by the present invention is described, such as Figure 4 shown.

[0082] Step 4: Obtain the dictionary vector corresponding to each iterative residual based on the rotation invariant subspace information.

[0083] Assume the initial state support set S 0 =φ, residual C 0 =R S , the sparse representation matrix

[0084] Let this be the kth iteration: the support set is S k , the residual is C k , the sparse representation matrix is The first P most significant solutions of the kth iteration are θ i ,i=1,2,...,P;

[0085] Let Span(A(θ)), Span(U S ) represent A(θ) and U S The space spanned by Span(A(θ))=Span(U S ), from the signal output equation R in step 3 S So Span(A(θ))=Span(R S )=Span(C 0 );

[0086] From the residual iteration equation, we know that Span(C k )=Span(R S );

[0087] From the above two equations, we can get Span(C k )=Span(A(θ)), so C k In the same subspace as A(θ);

[0088] From the above formula, we can get A(θ)·T=C k , where T is the rotation matrix.

[0089]

[0090] Among them A 1 (θ) and A 2 (θ) has N-1 rows; considering the number of signal sources P, A 1 (θ),A 2 (θ)∈C N-1×P , get A 2(θ) = A 1 (θ)·Ψ, where

[0091] Similarly, suppose

[0092]

[0093] From the above three formulas, we can get Available

[0094] From the above formula, we can get

[0095] From the above formula, we can get

[0096] Due to T -1 ΨT is similar to Ψ, so we can solve for Ψ;

[0097] The dictionary vector can be obtained:

[0098] Where i=1,2,...,P.

[0099] Step 5: Match the dictionary vectors and add them to the support set to calculate the updated residual. Let the index of the vector selected in the kth iteration be i k ,Pick will i k Add to S k , update the support set and sparse representation matrix:

[0100] S k+1 =S k ∪{i k},

[0101] The residual iteration equation of OMP (Orthogonal Matching Pursuit) is expressed as:

[0102]

[0103]

[0104] The subscript k represents the kth iteration, B represents the unknown number of the residual iteration equation, and B k represents the least squares solution of the kth iteration, C k represents the residual signal matrix of the kth iteration, R S represents the signal output equation obtained in step 3, represents the new dictionary matrix consisting of the k atoms selected at the kth iteration, express The pseudo-inverse matrix of . From the above formula, we can get the corresponding residual C k+1 , and iterate.

[0105] Step 6: Continue iterating until all the incoming wave directions are calculated. Repeat steps 4 and 5, and continue iterating until all the incoming wave directions are calculated.

[0106] It is to be understood that the present invention is described by some embodiments, and it is known to those skilled in the art that various changes or equivalent substitutions may be made to these features and embodiments without departing from the spirit and scope of the present invention. In addition, under the teachings of the present invention, these features and embodiments may be modified to adapt to specific circumstances and materials without departing from the spirit and scope of the present invention. Therefore, the present invention is not limited by the specific embodiments disclosed herein, and all embodiments falling within the scope of the claims of this application are within the scope of protection of the present invention.

Claims

1. A method for estimating a direction of arrival based on array signal processing, characterized in that: The following steps are involved: Step 1: Obtain incoming wave sampling data through the antenna array; Step 2, calculating the cyclic cross-correlation matrix by a cyclostationary algorithm; Step 3: Use singular value decomposition to obtain the corresponding signal subspace and construct a sparse reconstruction model; Step 4: Obtain the dictionary vector corresponding to each iteration residual based on the rotation invariant subspace information; Step 5: Match the dictionary vectors and add them to the support set to calculate the updated residual; Step 6: Continuously iterate until all incoming wave directions are calculated, repeat steps 4 and 5, and continuously iterate until all incoming wave directions are calculated; Step 1 specifically includes the following steps: using a uniform linear array to receive the signal, for a uniform linear array with N array elements spaced d apart, f0 is the signal frequency, λ0 is the signal wavelength; when M signal sources with a cyclic frequency of α are incident, where M≤N, the signal s received by the kth array element at the tth snapshot is obtained k (t), the signal wave direction θ i Satisfies -π / 2≤θ i ≤π / 2, where i=1,2,…,M, indicating that there are M signal sources in total; obtaining array sampled data S(t)=A(θ)X(t)+N(t), t=1,2,…,U, where A(θ) represents the received signal angle matrix, X(t) represents the basic signal matrix, N(t) represents the noise matrix, t represents the t-th snapshot data of the array, and it is assumed that a total of U snapshot data are received; The calculation method of step 2 includes the following steps: Step 2.1: Calculate the time delay τ between the kth array element and the i-th signal source ki ; Step 2.2: Based on the definition of cyclic cross-correlation and the array element receiving data, the array element s can be obtained. p and q Circular cross-correlation function of the received data Step 2.3: Where p,q=1,2,.,N, solve the cyclic cross-correlation function for N different array elements to obtain the cyclic cross-correlation matrix of the array received signal Where p,q represents the pth and qth array elements; Step 3 includes the following steps: For the circulant cross-correlation matrix obtained in step 2 Perform singular value decomposition to obtain unitary matrices U and V; arrange the singular values ​​in order from large to small, and divide the unitary matrices U and V into two parts: signal subspace and noise subspace according to the size of the singular values; span the signal subspace U from the eigenvectors of the M largest singular values S and V S ; Finally, the new signal output equation R is obtained S ; In step 3, the cyclic cross-correlation matrix The i-th column of is reconstructed like the sparse recovery model: Where i = 1, 2, ..., N; and r i Both said The i-th column of is the dictionary matrix, b i It is a dictionary Next i The sparse coefficient vector of ; represents the noise value, e i is an N×1 column vector, the i-th element is 1 and the rest are 0; according to Change it to matrix form as follows: Where B=[b1,b2,...,b N ] and b i (i=1,2,..,N) have the same sparse structure, that is, each non-zero element of the ideal value should appear in the same row of B, and the number of rows filled with 0 is the same, I N is the N×N identity matrix; in is an N×N matrix; The singular value decomposition of is expressed as: Where U and V are both unitary matrices, and Σ consists of singular values; therefore There are M large singular values; arrange the singular values ​​in order from large to small: λ1≥λ2≥...≥λ M ≥λ M+1 ≥...≥λ N ; The last NM singular values ​​are caused by noise and interference; according to the size of the singular values, the unitary matrices U and V are divided into two parts: where Σ S It is composed of the first M largest singular values, Σ N It is composed of the remaining singular values; Among them U S and V S represents the signal subspace spanned by the eigenvectors of the M largest singular values, and the Σ diagonal The remaining NM singular values ​​of constitute the noise subspace U N and V S ; A new signal output equation can be obtained: in is the noise interference error term; ΔR is given by and The error caused by the approximation of is a component of the signal subspace; Step 4 The calculation method includes the following steps: Step 4.1: Initialize the state support set S0 = φ, the residual C0 = R S , sparse representation matrix Step 4.2: Let this be the kth iteration: the support set is S k , the residual is C k , the sparse representation matrix is The first P most significant solutions of the kth iteration are θ i ,i=1,2,...,P; get A(θ)·T=C k , where T is the rotation matrix; let Where A(θ) represents the received signal angle matrix, A1(θ) represents the first N-1 rows of the matrix of A(θ), and A2(θ) represents the last N-1 rows of the matrix of A(θ). Represents C k The first N-1 rows of the matrix, Respectively represent C k The last N-1 rows of the matrix; Assume the transition matrix Ψ, where A2(θ) = A1(θ)·Ψ; calculate and get Due to T -1 ΨT is similar to Ψ and Ψ is a diagonal matrix, so Ψ can be solved; thus the dictionary vector α of the kth iteration is calculated i ,i=1,2,...,P; Step 5 includes the following steps: let the selected vector index of the kth iteration be i k ,Pick will i k Add to S k , update the support set and sparse representation matrix; thereby updating the residual iteration equation of orthogonal matching pursuit and obtaining a new residual solution C k+1 ; In step 5, the support set and sparse representation matrix are updated: S k+1 =S k ∪{i k }, The residual iteration equation of orthogonal matching pursuit is expressed as: The subscript k represents the kth iteration, B represents the unknown number of the residual iteration equation, and B k represents the least squares solution of the kth iteration, C k represents the residual signal matrix of the kth iteration, R S represents the signal output equation obtained in step 3, represents the new dictionary matrix consisting of the k atoms selected at the kth iteration, express The pseudo-inverse matrix, F refers to the two norm; from the above formula, the corresponding residual C k+1 , and iterate.

2. The method for estimating the direction of arrival based on array signal processing according to claim 1, characterized in that: In step 2.1, In step 2.2, in represents the cyclic autocorrelation of the basic signal of the ith signal source, τ pi represents the delay between the pth array element and the i-th signal source, τ qi represents the time delay between the qth array element and the i-th signal source, and for an equidistant uniform linear array, the following formula holds true: Since the noise cycle frequency is not α, its correlation coefficient at the α frequency is 0, and the cyclic cross-correlation function is: where s p (t) represents the data sampled from the pth array, s q (t) represents the data of the qth array sample; The cyclic cross-correlation matrix of the array received signal in step 2.3 is calculated as:

Images