CN104678350B - Arrival direction estimation based on TLS ESPRIT algorithms in extensive mimo system - Google Patents

Abstract

The present invention discloses a kind of method for realizing estimating two-dimensional direction-of-arrival using total least square ESPRIT algorithms in extensive mimo system, this method is first according to extensive mimo system receipt signal model, face battle array is divided into two partly overlapping submatrixs, it is met the array motion immovability required by ESPRIT algorithms；Then the signal subspace of subarray is calculated according to 1D ESPRIT algorithms, and then according to TLS criterions, estimates elevation information；Orientation angular estimation is finally realized using the direction matrix of a wherein submatrix, and is matched automatically with the elevation angle.The present invention efficiently solves arrival direction estimation problem in extensive mimo system using the thought of TLS ESPRIT algorithms, is not only able to realize azimuth and elevation angle pairing automatically, and computation complexity is relatively low.

Description

Arrival direction estimation based on TLS-ESPRIT algorithms in extensive mimo system

Technical field

The invention mainly relates to future mobile field.

Background technology

Following mobile data services amount will be exponentially increased, and traditional MIMO technology can not meet such demand. The end of the year 2010, AT＆T Labs scientist Thomas L.Marzetta propose extensive MIMO (Large Scale MIMO), Also it is Massive MIMO concept.Extensive MIMO is as a kind of key technology of future broadband wireless communication systems, in base station end The huge antenna of configured number, with depth digging utilization Spatial Dimension unlimited resources, lifting system spectrum efficiency and power.Big rule Mould mimo system is generally in base station end using face battle array, cylindrical array or special-shaped battle array.In extensive mimo system, due to base station Multiple antennas is configured, and causes base station acquisition channel condition information (Channel State Information, CSI) very difficult, is passed System base station precoding place one's entire reliance upon channel feedback no longer be applicable.2D (two dimension) DOA is estimated as channel precoding and provides essence True spatial information (si) so that the cost of implementation volume of optimal downlink precoding is greatly lowered.

In numerous extensive MIMO base station antenna configurations, URA is relatively easy, at the same can make full use of horizontal dimensions and The spatial information that vertical dimensions are provided, therefore, is more suitable for applying in extensive mimo system.Under URA configurations, letter is received Number simultaneously comprising vertical dimensions and horizontal dimensions information, 2D Mutual couplings seek to realize the orientation of incoming wave signal level dimension Angle and the elevation estimate of vertical dimension.

Mutual coupling is as one of important research content of array signal processing, in sonar, radar, radio communication system The fields such as system have important application value.At present, existing 1D DOA estimation theories have developed to obtain relative maturity, 1D DOA Method of estimation such as multiple signal classification (Multiple Signal Classification, MUSIC) algorithm, invariable rotary is empty Between algorithm (Estimating Signal Parameter via Rotational In variance Techniques, ESPRIT) algorithm is applied than wide, but is difficult to obtain practical application in two-dimensional estimation.Compared with 1D (one-dimensional) DOA, 2D DOA can preferably state signal space position, therefore research 2D DOA have more and are of practical significance.Utilize classical MUSIC algorithms , it is necessary to carry out spectrum peak search when carrying out two-dimensional estimation, operand is big, and consumes the substantial amounts of time.Two dimensional ESPRIT algorithm is preferable Ground solves this problem, but effective angle estimation can just be obtained by needing to match final result.On the one hand, by In antenna number increase and 2D DOA estimation in parameter dimension increase, cause the increase of estimation procedure complexity. If still using traditional subspace class algorithm such as 2D MUSIC, 2D ESPRIT algorithms estimation angle, it will face algorithm and answer The problems such as miscellaneous, computationally intensive and angle is matched, makes it cannot be used for actual direction-finding system；On the other hand, doing due to noise Disturb, sampled signal has error in itself.Existing most of DOA algorithm for estimating are directly carried out sampled signal as signal in itself Computing, does not meet actual direction finding environment.

The content of the invention

The deficiency that the present invention exists for current arrival direction estimation method, it is proposed that be based in extensive mimo system The arrival direction estimation method of TLS-ESPRIT algorithms, extensive MIMO is efficiently solved using the thought of TLS-ESPRIT algorithms Arrival direction estimation problem in system, is not only able to realize azimuth and elevation angle pairing automatically, and computation complexity is relatively It is low.

The present invention is to use URA in base station end, to receive different directions incoming wave.Believed according to the reception of extensive mimo system Number model, is divided into two partly overlapping submatrixs by face battle array first, it is met the array of ESPRIT algorithm requirements and moves constant Property；Then the signal subspace of subarray is calculated according to 1D ESPRIT algorithms；TLS criterions are recycled, in noise disturbance condition Under, estimate accurate elevation information；Orientation angular estimation is completed finally by the direction matrix of a wherein submatrix, and automatically with facing upward Match at angle.

The present invention comprises the following steps：

Step one：Base station uses M × N-dimensional URA, it is considered to which face battle array is located at XOZ planes, it is assumed that have K information source to incide this battle array Row, array received signal model is as follows：

X=AS+N

Wherein X is the dimension array received signal phasor of MN × 1, and S is the dimension space signal phasor of K × 1, and N is the Gauss white noise of MN × 1 Acoustic vector, A is that MN × K ties up direction matrix, and concrete form is as follows：

WhereinRepresent gram labor Roc inner product, a (u_i) be X-direction direction vector, a (v_i) be Z-direction direction vector.

Step 2：Consider that base station end aerial array is located at XOZ planes, prolong Z-direction and be divided into M × N-dimensional URA with shifting not Two submatrixs of denaturation, because the direction matrix of submatrix 1 and submatrix 2 differs a twiddle factor Φ, and the twiddle factor is only wrapped It is specific as follows containing elevation information：

Φ=diag [exp (- 2 π d cos θ₁/λ) exp(-2πd cosθ₂/λ) ... exp(-2πd cosθ_K/λ)]

Wherein d is array element spacing, and λ is velocity of wave, θ_i, i=1,2,3...K represent the elevation angle of i-th of user.

Step 3: according to the signal subspace E of 1D ESPRIT algorithm computing arrays_s, define selection matrix J₁,J₂, therefrom Select the corresponding signal subspace of two submatrixsDetailed process is as follows：

Definition

Then have：

Then submatrix 1 and the corresponding signal subspace of submatrix 2 are：

Step 4: considering there is noise disturbance, according to TLS criterions, elevation information is estimated, is specifically expressed as follows：

According to TLS criterions, definition:

It is rightCarry out Eigenvalues Decomposition：

And E is resolved into K × K submatrix, i.e.,：

CalculateAnd Eigenvalues Decomposition is carried out to Ψ：

Ψ=H Ω H^-1=T Φ T^-1

Eigenvalue λ in Ω_i, i=1, the diagonal entry in 2...K homographies Φ, H=T calculates elevation estimate value：

θ_i=arccos { angle (λ_i)λ/2πd}

Wherein d is array element spacing, and λ is velocity of wave, θ_i, i=1,2,3...K represent the elevation angle of i-th of user, ()^HExpression is asked Take the conjugate transposition of matrix, ()^-1Inverse of a matrix is asked in expression.

Step 5: the direction matrix of submatrix 1 can be by signal subspaceWith non-singular matrix H product representation, its side is utilized To matrix column vector can estimation orientation angle, specific formula is as follows：

In i-th, i=1,2...K column vectorsMeet：

WhereinRepresent that i-th arranges preceding 1:M (N-2) individual element,Represent N-1 after the i-th row:M(N-1) Individual element.

OrderAzimuth φ can be estimated according to above formula：

φ_i=arccos (angle (b_i)λ/(2πd sinθ_i))

The present invention considers the noise jamming of sampled signal presence, using TLS criterions, estimates the accurate elevation angle；Pass through side Orientation angular estimation is completed to rectangular array Vector rotation invariant relation.The present invention realizes that 2D angles are estimated merely with 1D DOA algorithm for estimating Meter, significantly reduces computation complexity, and can realize angle automatic matching.

Brief description of the drawings

Fig. 1 is the FB(flow block) of 2D DOA methods in extensive mimo system proposed by the present invention；

Fig. 2 is extensive MIMO 3D model schematics；

Fig. 3 is that the URA that base station end is used divides schematic diagram；

Fig. 4 is the arrival direction estimation method based on TLS-ESPRIT algorithms in extensive mimo system proposed by the present invention Flow chart.

Embodiment

In order that the object, technical solutions and advantages of the present invention are clearer, below in conjunction with accompanying drawing the present invention is made into The detailed description of one step：

Fig. 1 show the FB(flow block) of 2D DOA methods in extensive mimo system proposed by the present invention, base station end configuration M × N-dimensional URA, wherein M represent that URA levels tie up antenna number, and N represents URA vertical dimension antenna numbers, it is ensured that each transmitting antenna can Enough to handle horizontal dimensions and the information in vertical dimensions simultaneously, it is K that user terminal, which sends signal,.According to extensive mimo system Receipt signal model, is divided into two partly overlapping submatrixs by face battle array first, subarray is calculated by 1D ESPRIT algorithms Signal subspace；Then TLS criterions are utilized, elevation information is estimated；Finally using the direction matrix of a wherein submatrix come the side of completion Parallactic angle is estimated, and realizes automatic matching.

Fig. 2 is extensive MIMO 3D model schematics, and base station is highly h, and user is highly h_u, base station end, which is configured, to be located at The M row N row URA of X-Z plane, array element spacing is d, and it is K to send signal number.There is an orientation for each transmission signal Angle φ and elevation angle theta are corresponded to therewith, and system model is as follows：

X=AS+N

Wherein X is the dimension array received signal phasor of MN × 1, and S is the dimension space signal phasor of K × 1, and N is the Gauss white noise of MN × 1 Acoustic vector, A is that MN × K is tieed up under direction matrix, concrete form：

a(u_i)=[1 exp (- j2 π d sin θs_icosφ_i/λ) ... exp(-j(M-1)2πd sinθ_icosφ_i/λ)]^T

a(v_i)=[1 exp (- j2 π d cos θ_i/λ) ... exp(-j(N-1)2πd cosθ_i/λ)]^T

WhereinRepresent gram labor Roc inner product, a (u_i) be X-direction direction vector, a (v_i) be Z-direction direction vector, θ_iAnd φ_iThe elevation angle and azimuth of i-th of user is represented respectively, and i=1,2,3...K, d is array element spacing, and λ is velocity of wave, ()^TTable Show the transposition for asking for matrix.

Fig. 3 is that the URA that base station end is used divides schematic diagram, it is considered to which base station end aerial array is located at XOZ planes, prolongs Z-direction Planar array is divided into two submatrixs with motion immovability, two submatrix array numbers are identical, the direction square of submatrix 1 and submatrix 2 Battle array one twiddle factor Φ of difference, the twiddle factor only includes elevation information θ, and concrete form is as follows：

A₁Φ=A₂

Φ=diag [exp (- 2 π d cos θ₁/λ) exp(-2πd cosθ₂/λ) ... exp(-2πd cosθ_K/λ)]

Wherein A₁,Α₂The respectively direction matrix of submatrix 1 and submatrix 2, d is array element spacing, and λ is velocity of wave, θ_i, i=1,2, 3...K the elevation angle of i-th of user is represented.

Fig. 4 is the arrival direction estimation method flow diagram based on TLS-ESPRIT algorithms in extensive mimo system, specific real Existing step is as follows：

Step 41, as shown in figure 3, two partly overlapping submatrixs, submatrix 1 will be divided into along Z axis positioned at the URA in X-Z faces A twiddle factor Φ is differed with the direction matrix of submatrix 2, the twiddle factor only includes elevation information θ.

Step 42, according to 1D ESPRIT algorithms, the covariance matrix that computing array is exported first：

Wherein N is the fast umber of beats of signal sampling, then to R_XXEigenvalues Decomposition is carried out, signal subspace E is obtained_s：

Wherein σ²For noise variance, ()^HThe conjugate transposition of matrix is asked in expression.

Step 43, selection matrix J is defined₁,J₂Signal subspace E_sIn select the corresponding signal subspace of two submatrixsSpecifically it is expressed as follows：

Definition：

Then have：

WhereinRepresent gram labor Roc inner product.So corresponding signal subspace of submatrix is：

Step 44, the elevation angle is estimated according to TLS criterions, detailed process is as follows：

Under A, ideal conditions, the signal subspace E of submatrix 1 and submatrix 2₁,E₂With corresponding direction matrix A₁,Α₂Zhang Chengxiang As subspace, then the non-singular matrix T of existence anduniquess so that following formula set up：

A₁=E₁T A₂=E₂T

According to conclusion A in Fig. 3₁Φ=A₂, E can be obtained₁TΦT^-1=E₂；Define Ψ=T Φ T^-1, then have E₁Ψ=E₂.Wherein, ()^-1Inverse of a matrix is asked in expression.

There is error in B, actual samples signal, defence the signal subspace that difference is solved jointlyAll there is noise disturbance：

Ψ is solved using TLS criterions more more accurate than traditional ESPRTI algorithms, be specifically expressed as follows：

Definition：

It is rightCarry out Eigenvalues Decomposition：

And E is resolved into K × K submatrix, i.e.,

CalculateAnd Eigenvalues Decomposition is carried out to Ψ：

Ψ=H Ω H^-1=T Φ T^-1

Eigenvalue λ in Ω_i, i=1, the diagonal entry in 2...K homographies Φ, H=T calculates elevation estimate value：

θ_i=arccos { angle (λ_i)λ/2πd}

Wherein d is array element spacing, and λ is velocity of wave, ()^HThe conjugate transposition of matrix, () are asked in expression^-1Matrix is asked in expression It is inverse.

Step 45, using the direction matrix estimation orientation angle of submatrix 1, the corresponding direction matrix of submatrix 1 can be empty by signal subspace Between product representation with above-mentioned non-singular matrix H, be specifically expressed as follows：

The direction matrix A of submatrix 1₁：

A₁=[α₁ α₂ … α_i … α_K]_M(N-1)×K

Wherein α_iConcrete form be：

By α_iConcrete form can be seen that A₁In i-th arrange preceding 1:M (N-2) individual element and rear N-1:M (N-1) individual element is expired Sufficient invariable rotary relation：

Utilize the signal subspace of submatrix 1Go out the corresponding direction matrix of submatrix with non-singular matrix H product estimation：

Step 46, utilizeThe invariable rotary relation of middle column vector is by twiddle factorSolve, estimation orientation Angle φ, specific formula is as follows：

φ_i=arccos (angle (b_i)λ/(2πd sinθ_i))

Wherein d is array element spacing, and λ is velocity of wave, φ_i, i=1,2,3...K represent the azimuth of i-th of user, ()⁺Represent Ask for the pseudoinverse of matrix.

Above is the preferred embodiments of the present invention are described in detail, it will be appreciated that preferred embodiment is only for saying The bright present invention, the protection domain being not intended to be limiting of the invention.

Claims (3)

1. the arrival direction estimation method based on TLS-ESPRIT algorithms in extensive mimo system, it is characterised in that including step Suddenly：Base station end is using M × N-dimensional uniform rectangular face battle array URA, according to the receipt signal model of extensive mimo system, first by face Battle array is divided into two partly overlapping submatrixs, the array motion immovability for making it meet ESPRIT algorithm requirements；And then according to one-dimensional ESPRIT algorithms calculate the signal subspace of submatrix；Then TLS criterions are utilized, elevation estimate is realized；Finally utilize wherein one son The direction matrix of battle array estimates corresponding azimuth；

Described two partly overlapping submatrixs are submatrix 1 and submatrix 2, two submatrixs with motion immovability, two submatrix array numbers Identical, the direction matrix of submatrix 1 and submatrix 2 differs a twiddle factor Φ, and the twiddle factor only includes elevation information θ, specifically Form is as follows：

A₁Φ=A₂

Φ=diag [exp (- 2 π dcos θ₁/λ) exp(-2πdcosθ₂/λ) ... exp(-2πdcosθ_K/λ)]

Wherein A₁,Α₂The respectively direction matrix of submatrix 1 and submatrix 2, d is array element spacing, and λ is velocity of wave, θ_i, i=1,2,3...K Represent the elevation angle of i-th of user；

It is described according to 1D ESPRIT algorithms calculate submatrix signal subspace be：According to one-dimensional ESPRIT algorithms, submatrix is sampled The covariance matrix of data carries out Eigenvalues Decomposition, obtains signal subspace E_s, define selection matrix J₁,J₂, from E_sIn select with Submatrix 1 and the corresponding signal subspace of submatrix 2WithI.e.：

2. arrival direction estimation method according to claim 1, it is characterised in that the utilization TLS criterions, realizes the elevation angle Estimation is：Consider under conditions of noise disturbance, according to TLS criterions, estimate accurate elevation information, definitionIt is rightCarry out Eigenvalues Decomposition and obtain eigenmatrix E, and E is resolved into K × K Matrix, i.e.,：

$E=\left[\begin{array}{cc}{E}_{11}& E_{12}\\ E_{21}& E_{22}\end{array}\right]$

CalculateAnd Eigenvalues Decomposition is carried out to Ψ：

Ψ=H Ω H^-1=T Φ T^-1

Eigenvalue λ in Ω_i, i=1, the diagonal entry in 2...K homographies Φ, H=T calculates elevation estimate value：

θ_i=arccos { angle (λ_i)λ/2πd}

Wherein d is array element spacing, and λ is velocity of wave, θ_iRepresent the elevation angle of i-th of signal, ()^HThe conjugate transposition of matrix is asked in expression, ()^-1Inverse of a matrix is asked in expression.

3. arrival direction estimation method according to claim 1, it is characterised in that the direction using a wherein submatrix The corresponding azimuth specific method of Matrix Estimation is：The direction matrix of submatrix 1 can be by its signal subspaceWith non-singular matrix H's Product representation：The invariable rotary relation of utilization orientation matrix column vector can estimate corresponding azimuth φ：

φ_i=arccos (angle (b_i)λ/(2πdsinθ_i))

Wherein RepresentI-th arranges preceding 1:M (N-2) individual element,Represent N-1 after the i-th row:M (N-1) individual element, d is array element spacing, and λ is velocity of wave, θ_iRepresent facing upward for i-th signal Angle, φ_iRepresent the azimuth of i-th of signal, ()⁺The pseudoinverse of matrix is asked in expression.