CN108710102B - Direction-of-arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of co-prime array - Google Patents

Direction-of-arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of co-prime array Download PDF

Info

Publication number
CN108710102B
CN108710102B CN201810460228.7A CN201810460228A CN108710102B CN 108710102 B CN108710102 B CN 108710102B CN 201810460228 A CN201810460228 A CN 201810460228A CN 108710102 B CN108710102 B CN 108710102B
Authority
CN
China
Prior art keywords
signal
array
virtual
angle
fourier transform
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201810460228.7A
Other languages
Chinese (zh)
Other versions
CN108710102A (en
Inventor
张宗煜
史治国
周成伟
陈积明
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University ZJU
Original Assignee
Zhejiang University ZJU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang University ZJU filed Critical Zhejiang University ZJU
Priority to CN201810460228.7A priority Critical patent/CN108710102B/en
Publication of CN108710102A publication Critical patent/CN108710102A/en
Application granted granted Critical
Publication of CN108710102B publication Critical patent/CN108710102B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/02—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using radio waves
    • G01S3/14—Systems for determining direction or deviation from predetermined direction
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/78—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using electromagnetic waves other than radio waves
    • G01S3/782—Systems for determining direction or deviation from predetermined direction
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/80—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves
    • G01S3/802—Systems for determining direction or deviation from predetermined direction
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
    • G06F17/141—Discrete Fourier transforms
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • General Engineering & Computer Science (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Algebra (AREA)
  • Discrete Mathematics (AREA)
  • Computing Systems (AREA)
  • Electromagnetism (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)
  • Radar Systems Or Details Thereof (AREA)
  • Variable-Direction Aerials And Aerial Arrays (AREA)

Abstract

The invention discloses a direction of arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of a co-prime array, which mainly solves the problems that the existing method is high in calculation complexity and cannot estimate signal power at the same time. The method comprises the following implementation steps: a receiving end constructs a co-prime array; receiving an incident signal by utilizing a co-prime array and modeling; deriving a second order equivalent virtual signal corresponding to the augmented virtual uniform linear array from the co-prime array received signal; defining an angle-space frequency and describing a second-order equivalent virtual signal of the virtual uniform linear array by the angle-space frequency; performing inverse discrete Fourier transform on the second-order equivalent virtual signal described by the angle-space frequency to construct a space power spectrum; and performing spectral peak search according to the constructed spatial power spectrum to obtain the direction of arrival estimation and power estimation information of the signal. The invention improves the performance of the degree of freedom of the direction of arrival estimation, reduces the computational complexity of the direction of arrival estimation and can simultaneously obtain the power estimation information of the signal.

Description

Direction-of-arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of co-prime array
Technical Field
The invention belongs to the technical field of signal processing, particularly relates to statistical signal processing of radar signals, acoustic signals and electromagnetic signals, and particularly relates to a direction of arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of a co-prime array, which can be used for passive positioning and target detection.
Background
Direction-of-Arrival (DOA) estimation is a basic problem in the field of array signal processing, and it is a method of receiving signals by using a sensor array and statistically processing the received signals by a series of signal processing methods to obtain Direction-of-Arrival information contained in the signals, and is widely applied in the fields of radar, sonar, voice, wireless communication, and the like.
The degree of freedom of the DOA estimation method refers to the number of incident signal sources that it can estimate. The traditional DOA estimation method generally adopts a uniform linear array for signal receiving and modeling, but the degree of freedom of the DOA estimation method based on the uniform linear array is limited by the number of physical array elements. When the number of incident signal sources is larger than the number of physical array elements in the array, the DOA estimation method based on the uniform linear array cannot obtain an effective estimation result. The co-prime array is a non-uniform sparse array with a systematic structure, can break through the limitation of the degree of freedom of the traditional uniform linear array, and realizes the improvement of the degree of freedom performance under the condition of a certain number of physical array elements, so that the co-prime array has been widely paid attention to in recent years in the academic world. The principle of the DOA estimation method based on the co-prime array is that the co-prime array is deduced to a virtual domain by utilizing the property of prime numbers, and a second-order equivalent virtual signal corresponding to a virtual uniform linear array is constructed for DOA estimation. The number of virtual array elements in the virtual uniform linear array is larger than that of physical array elements, so that the degree of freedom is improved.
Most of existing DOA estimation methods based on the co-prime array perform complex operations including complex matrix operations such as inversion and eigenvalue decomposition and processes such as design and solution of convex optimization problems on the basis of the second-order equivalent virtual signals of the co-prime array, the operation processes result in higher calculation complexity, certain challenges are faced in application scenes with higher real-time requirements, and hardware in actual systems is difficult to implement. In addition, many existing DOA estimation methods cannot obtain power information of a signal at the same time as obtaining estimation information of a direction of arrival. However, the signal power is also an important parameter for describing the signal, and is of great significance for detecting and identifying the target.
Disclosure of Invention
The invention aims to provide a direction of arrival estimation method based on the inverse discrete Fourier transform of a second-order equivalent virtual signal of a co-prime array, aiming at the defects in the prior art, the DOA estimation freedom performance is improved by performing the inverse discrete Fourier transform on the second-order equivalent virtual signal of the co-prime array described by adopting angle-space frequency, and the direction of arrival estimation of the signal and corresponding signal power information are obtained at the same time. The method provided by the invention has lower computational complexity and is easy to realize in hardware in an actual system.
The purpose of the invention is realized by the following technical scheme: a direction of arrival estimation method based on inverse discrete Fourier transform of a second-order equivalent virtual signal of a co-prime array comprises the following steps:
(1) the receiving end uses 2M + N-1 array elements and constructs according to a co-prime array structure; wherein M and N are relatively prime integers;
(2) suppose there are L from [ theta ═ theta1,θ2,…,θL]TDirectional far-field narrow-band incoherent signal source [. ]]TRepresenting a transposition operation, receiving an incident signal by using the co-prime array constructed in step (1), and modeling a co-prime array received signal x (t) at time t as follows:
Figure BDA0001660501050000021
wherein x (t) is a (2M + N-1) × 1-dimensional vector, sl(t) is the waveform of the first incident signal, n (t) is the noise component independent of each signal source, a (theta)l) To correspond to thetalThe guiding vector of the co-prime array of the directional signal source can be expressed as
Figure BDA0001660501050000022
Wherein, mu i1,2,3, …,2M + N-1 represents the actual position of the ith physical array element in the co-prime array, and the position of the first physical array element is mu1λ is the wavelength of the incident narrowband signal, and j is the imaginary unit.
Constructing a covariance matrix according to the co-prime array received signal x (t): sampling covariance matrix obtained using T sampling snapshots
Figure BDA0001660501050000023
For the theoretical covariance matrix RxCarrying out approximate substitution;
(3) deriving a second order equivalent virtual signal corresponding to the augmented virtual uniform linear array from the co-prime array received signal: sampling covariance matrix by vectorization
Figure BDA0001660501050000031
Obtaining a virtual array equivalent virtual signal y:
Figure BDA0001660501050000032
wherein,
Figure BDA0001660501050000033
is (2M + N-1)2The steering matrix is maintained at × L,
Figure BDA0001660501050000034
is a vector containing the L incident signal source powers,
Figure BDA0001660501050000035
for noise power, I ═ vec (I)2M+N-1)，I2M+N-1Representation (2M + N-1 × (2M + N-1 dimensional identity matrix, vec (-) represents vectorization operations, i.e., stacking columns in the matrix in order to form a new vector (-)*It is meant a conjugate operation of the two,
Figure BDA0001660501050000036
represents the kronecker product;
non-uniform virtual array S corresponding to vector yDExpressed as:
SD＝{±(Mn-Nm)d,0≤n≤N-1,0≤m≤2M-1}，
where d is half the wavelength of the incident narrowband signal, i.e.
Figure BDA0001660501050000037
Selecting a non-uniform virtual array SDThe virtual array elements of the medium-maximum continuous part form a virtual uniform linear array S containing 2V +1 virtual array elementsV＝{-Vd,-(V-1)d,…,0,…,(V-1)d,Vd}，V＝MN+M-1；
Selecting the vector y corresponding to SVThe equivalent virtual signals of the positions of the virtual array elements form a second-order equivalent virtual signal z corresponding to the virtual uniform linear arrayθIt can be expressed as:
Figure BDA0001660501050000038
wherein B (θ) ═ B (θ)1),b(θ2),…b(θL)]Column l thereof
Figure BDA0001660501050000039
To correspond to thetalA virtual uniform linear array steering vector of the directional signal source, e is selected from i corresponding to SVVectors composed of elements of the medium array elements;
(4) the angle-spatial frequency is defined and used to describe the second order equivalent virtual signal of the virtual uniform line array. The angle-space frequency is defined as the number of signal cycles of a narrow-band signal from the direction theta in space, which is propagated within the propagation distance difference between adjacent array elements. Second-order equivalent virtual signal z of the virtual uniform linear array in the step (3)θCan be equivalently expressed in the angular-spatial frequency domain as:
Figure BDA0001660501050000041
wherein, B (ξ) ═ B (ξ)1),b(ξ2),…b(ξL)]Column l thereof
Figure BDA0001660501050000042
Figure BDA0001660501050000043
ξ＝[ξ1,ξ2,…,ξL]TAngle-space frequencies corresponding to L angles contained in θ;
(5) performing inverse discrete Fourier transform on the second-order equivalent virtual signal described by the angle-space frequency, and constructing a space power spectrum: second order equivalent virtual signal z represented by angle-space frequency by inverse discrete Fourier transformξConverting the space domain to obtain a K × 1 dimensional space response psi;
constructing a spatial power spectrum, wherein the horizontal axis of the spectrum represents an angle theta, and the relation between the angle theta and the k-th element of the spatial response vector can be expressed as:
Figure BDA0001660501050000044
wherein K is 0,1, …, K-1, arccos (·) is an inverse cosine function, and h is a guarantee
Figure BDA0001660501050000045
Coefficients satisfying the domain of the inverse cosine function when
Figure BDA0001660501050000046
When h is-1, when
Figure BDA0001660501050000047
When h is 0; the vertical axis of the spectrum represents the modulus p (k) of the kth element in the spatial response vector;
(6) and estimating the direction of arrival and the signal power according to the obtained spatial power spectrum. And (5) performing spectrum peak searching operation on the spatial power spectrum in the step (5), wherein the angle corresponding to the first L peak values with the maximum amplitude is the direction of arrival estimation of the L incident signals, and the peak value amplitude is the power estimation value of the corresponding signal.
Further, the relatively prime array structure in step (1) can be specifically described as follows: first, a set of coprime integers M, N is selected to construct a pair of sparse uniform linear sub-arrays. The first subarray comprises 2M Nd-spaced array elements at positions 0, Nd, …, (2M-1) Nd; the second subarray comprises N array elements spaced Md apart at positions 0, Md, …, (N-1) Md. And then combining the two sub-arrays according to a mode that the first array elements are overlapped to obtain a non-uniform co-prime array containing 2M + N-1 physical array elements.
Further, in step (3), if there are a plurality of different elements in y and SVIf the same virtual array element position corresponds to each other, any one of the elements is selected as a composition vector zθOf (2) is used.
Further, the angle-space frequency in step (4) is a frequency related to the incident angle of the signal, and is specifically defined by the following way: the arrival of a narrow-band signal from the direction θ at an adjacent array element with a distance d in space generates a propagation distance difference u, which can be expressed as:
u＝dcosθ。
in the case of a fixed array element pitch, the propagation distance difference u varies only with the incident signal angle θ, and therefore the definition of the angle-spatial frequency in step (4) is obtained. The angle-spatial frequency ξ versus incident signal angle θ can be expressed as:
Figure BDA0001660501050000051
further, the second-order equivalent virtual signal z obtained in the step (5) through inverse discrete fourier transformξThe K × 1 dimensional spatial response ψ of can be expressed as:
Figure BDA0001660501050000052
wherein,
Figure BDA0001660501050000053
denotes the inverse discrete Fourier transform operation, FKIs a K-point inverse discrete fourier transform matrix, which can be expressed as:
Figure BDA0001660501050000054
the resulting spatial response ψ is a K × 1 dimensional vector.
Further, in the step (5), the constructed spatial spectrum reflects the response amplitude at each angle in space, where there are L peaks corresponding to L incident signals.
Further, there are L peaks in the spatial power spectrum corresponding to the L incident signals in the step (5), and the conclusion is obtained as follows: the relation p (k) used to construct the spatial power spectrum may be specifically expressed as:
Figure BDA0001660501050000055
wherein,
Figure BDA0001660501050000056
for the power of the l-th signal, (-) represents an impulse function, which is an integer, representing the impulse function as a sequence of periodic impulse strings. According to the nature of the impulse function, only when
Figure BDA0001660501050000057
Figure BDA0001660501050000058
When P (k) is the peak value. Since the L incident signals have different angle-space frequencies, each of which can and can only cause one peak in the spatial power spectrum, there are L peaks in the spatial power spectrum corresponding to the L incident signals.
Compared with the prior art, the invention has the following advantages:
(1) the method provided by the invention obtains the space response by carrying out inverse discrete Fourier transform on the second-order equivalent virtual signal of the cross-prime array described by the angle-space frequency, constructs the space power spectrum based on the space response, obtains the estimation result of the direction of arrival by searching the spectral peak of the constructed space power spectrum, avoids the complex calculation processes of design solution of a convex optimization problem, matrix inversion, matrix eigenvalue decomposition and the like commonly used in the traditional direction of arrival estimation method, and can simultaneously obtain the power estimation of the signal on the basis of ensuring the improvement of DOA estimation freedom degree performance.
(2) The method provided by the invention adopts the inverse discrete Fourier transform to estimate the direction of arrival, effectively reduces the computational complexity, better meets the application requirement with higher real-time requirement, and is easy to realize in hardware in an actual system.
Drawings
FIG. 1 is a block diagram of the overall flow of the method of the present invention;
FIG. 2 is a schematic diagram of a pair of sparse uniform subarrays constituting a co-prime array according to the present invention;
FIG. 3 is a schematic diagram of the structure of the co-prime array of the present invention;
FIG. 4 is a schematic diagram of the signal power estimation accuracy at different SNR for the proposed method embodying the present invention;
FIG. 5 is a schematic diagram of signal power estimation accuracy at different sampling fast-beat numbers for the proposed method embodying the present invention;
fig. 6 is a schematic diagram of a spatial power spectrum for embodying the degree of freedom performance of the proposed method.
Detailed Description
The technical means and effects of the present invention will be described in further detail below with reference to the accompanying drawings.
The relatively prime array can perform statistical processing on the second-order equivalent virtual array signals, so that the degree of freedom performance is improved under the condition of a certain physical array element number, and the relatively prime array has attracted wide attention in recent years in academic circles. In the application of an actual system, the existing DOA estimation method has high calculation complexity, is difficult to meet the application scene with high real-time requirement, and has certain difficulty in realizing the hardware in the actual system in the complex calculation process. In addition, most existing DOA estimation methods cannot obtain the direction of arrival estimation and estimate the power of each signal source at the same time. In view of the above problems, the present invention provides a direction of arrival estimation method based on the inverse discrete fourier transform of a second-order equivalent virtual signal of a co-prime array, and referring to fig. 1, the implementation steps of the present invention are as follows:
the method comprises the following steps: and 2M + N-1 actual array elements are used at a receiving end to construct a co-prime array. First, a set of coprime integers M, N is selected to construct a pair of sparse uniform linear sub-arrays. The first subarray comprises 2M Nd-spaced array elements at positions 0, Nd, …, (2M-1) Nd; the second sub-array comprises N array elements with a spacing Md and the positions are 0, Md, …, (N-1) Md, wherein the unit spacing d is the half wavelength of the incident narrowband signal, i.e. d ═ λ/2. And then combining the two sub-arrays according to a mode that the first array elements are overlapped to obtain a non-uniform co-prime array containing 2M + N-1 actual array elements.
Step two: a co-prime array is used to receive and model the signal. Suppose there are L from [ theta ═ theta1,θ2,…,θL]TDirectional far-field narrow-band incoherent signal source [. ]]TRepresenting a transposition operation, receiving an incident signal by using the relatively prime array constructed in the step one, and modeling a relatively prime array received signal x (t) at a time t as follows:
Figure BDA0001660501050000071
wherein x (t) is a (2M + N-1) × 1-dimensional vector, sl(t) is the waveform of the first incident signal, n (t) is the noise component independent of each signal source, a (theta)l) To correspond to thetalThe guiding vector of the co-prime array of the directional signal source can be expressed as
Figure BDA0001660501050000072
Wherein, mu i1,2,3, 2M + N-1 tableShowing the actual position of the ith physical array element in the co-prime array, and the position of the first physical array element is mu1λ is the wavelength of the incident narrowband signal, and j is the imaginary unit. In practice, a sampling covariance matrix obtained by using T sampling snapshots
Figure BDA0001660501050000073
For the theoretical covariance matrix RxThe approximate substitution is made and,
Figure BDA0001660501050000074
and RxCan be respectively represented as
Figure BDA0001660501050000075
Rx＝E[x(t)xH(t)]，
Wherein, (.)HDenotes the conjugate transposition, E [. cndot.)]Representing a mathematical expectation.
Step three: a second order equivalent virtual signal corresponding to the augmented pseudo-uniform linear array is derived from the co-prime array received signal. The virtual array equivalent virtual signal y can be quantized by the sampling covariance matrix in the second step
Figure BDA0001660501050000081
Obtaining:
Figure BDA0001660501050000082
wherein,
Figure BDA0001660501050000083
is (2M + N-1)2The steering matrix is maintained at × L,
Figure BDA0001660501050000084
is a vector containing the L incident signal source powers,
Figure BDA0001660501050000085
for noise power, I ═ vec (I)2M+N-1)，I2M+N-1Representing a (2M + N-1) × (2M + N-1) dimensional identity matrix, vec (·) representing a vectorization operation, i.e., stacking columns in the matrix in sequence to form a new vector, (·)*It is meant a conjugate operation of the two,
Figure BDA0001660501050000086
representing the kronecker product. Non-uniform virtual array S corresponding to vector yDCan be expressed as:
SD＝{±(Mn-Nm)d,0≤n≤N-1,0≤m≤2M-1}，
where d is half the wavelength of the incident narrowband signal, i.e.
Figure BDA0001660501050000087
Selecting a non-uniform virtual array SDThe virtual array elements of the medium-maximum continuous part form a virtual uniform linear array S containing 2V +1 virtual array elementsV{ -Vd, - (V-1) d, …,0, …, (V-1) d, Vd }, wherein V ═ MN + M-1. Selecting the vector y corresponding to SVThe equivalent virtual signals of the positions of the virtual array elements form a second-order equivalent virtual signal z corresponding to the virtual uniform linear arrayθIt can be expressed as:
Figure BDA0001660501050000088
wherein B (θ) ═ B (θ)1),b(θ2),…b(θL)]Column l thereof
Figure BDA0001660501050000089
To correspond to thetalA virtual uniform linear array steering vector of the directional signal source, e is selected from i corresponding to SVThe elements of the array elements constitute a vector. If there are a plurality of different elements in y and SVIf the same virtual array element position corresponds to each other, any one of the elements is selected as a composition vector zθOf (2) is used.
Step four: the angle-spatial frequency is defined and used to describe the second order equivalent virtual signal of the virtual uniform line array. The angle-space frequency is defined as the number of signal cycles of a narrow-band signal from the direction theta in space, which is propagated within the propagation distance difference between adjacent array elements. Specifically, a narrow-band signal from the direction θ in space arrives at an adjacent array element with a distance d, which generates a propagation distance difference, which can be expressed as:
u＝dcosθ,
under the condition that array elements are spaced at a certain distance, the propagation distance difference varies with the incident angle theta, and the angle-space frequency is defined as the number of periods of narrowband signal propagation on the propagation distance difference, namely:
Figure BDA0001660501050000091
correspondingly, the second-order equivalent virtual signal z of the virtual uniform linear array in the step threeθCan be equivalently expressed in the angular-spatial frequency domain as:
Figure BDA0001660501050000092
wherein ξ ═ ξ1,ξ2,…,ξL]TAngle-spatial frequencies corresponding to the L angles contained in θ.
Step five: and performing inverse discrete Fourier transform on the second-order equivalent virtual signal described by the angle-space frequency to obtain a space response. The second-order equivalent virtual signal z represented by the angle-spatial frequency can be represented by an inverse discrete fourier transformξThe transformation into the spatial domain, and hence its spatial response ψ, can be expressed as:
Figure BDA0001660501050000093
wherein,
Figure BDA0001660501050000094
denotes the inverse discrete Fourier transform operation, FKIs a K-point inverse discrete fourier transform matrix, which can be expressed as:
Figure BDA0001660501050000095
the resulting spatial response ψ is a K × 1 dimensional vector. Constructing a spatial power spectrum, wherein the horizontal axis of the spectrum represents an angle theta, and the relation between the angle theta and the k-th element of the spatial response vector can be expressed as:
Figure BDA0001660501050000096
wherein K is 0,1, …, K-1, arccos (·) is an inverse cosine function, and h is a guarantee
Figure BDA0001660501050000097
Coefficients satisfying the domain of the inverse cosine function when
Figure BDA0001660501050000098
When h is-1, when
Figure BDA0001660501050000099
When h is 0; the vertical axis of the spectrum represents the modulus p (k) of the kth element in the spatial response vector, which can be expressed as:
P(k)＝|[ψ]k|，
wherein [ ·]kRepresenting the kth element in the vector, |, represents the modulus of the complex number. Specifically, p (k) can be represented as:
Figure BDA0001660501050000101
wherein,
Figure BDA0001660501050000102
for the power of the l-th signal, (-) represents an impulse function, which is an integer, representing the impulse function as a sequence of periodic impulse strings. According to the nature of the impulse function, only when
Figure BDA0001660501050000103
When P (k) is the peak value. Due to L piecesThe incident signals have different angle-space frequencies, each of which can and cannot cause only one peak in the spatial power spectrum, and thus there are L peaks in the spatial power spectrum corresponding to the L incident signals.
Step six: and estimating the direction of arrival and the signal power according to the obtained spatial power spectrum. And D, performing spectrum peak searching operation on the spatial power spectrum in the step five, and arranging the peak values of the spatial power spectrum from high to low, wherein the maximum first L peak values correspond to L incident signals, and the peak values of the L spectrum peaks are power estimated values of corresponding signals respectively.
The direction of arrival estimation method provided by the invention can simultaneously acquire direction of arrival estimation information and power estimation information of corresponding signals through a spectrum peak search process of the space power spectrum. Compared with the traditional direction of arrival estimation method based on the uniform linear array, the method provided by the invention ensures that the performance of the degree of freedom of the direction of arrival estimation is improved, and the calculation complexity is only
Figure BDA0001660501050000104
The method well meets the actual application scene with higher requirements on real-time performance, and is easier to realize in hardware in an actual system.
The effect of the present invention will be further described with reference to the simulation example.
Simulation example 1: the incident signal is received by a co-prime array, and the parameters are selected to be M-9 and N-10, that is, the co-prime array comprises 2M + N-1-27 physical array elements. The fixed sampling fast beat number T is 1000, and under the condition of respectively testing different signal-to-noise ratios, the power estimation accuracy of a single signal source and 30 signal sources is expressed by Root Mean Square Error (RMSE), and the calculation formula is as follows:
Figure BDA0001660501050000105
whereinQ is the number of Monte Carlo trials,
Figure BDA0001660501050000111
represents the power estimate of the ith signal in the qth monte carlo test,
Figure BDA0001660501050000112
is the actual power of the l-th signal. In this simulation, Q is 1000. For the case of a single signal source, the incident direction of the signal source in each Monte Carlo test satisfies the Gaussian distribution
Figure BDA0001660501050000113
In the case of 30 signal sources, the incident directions of the signal sources are uniformly distributed in the spatial angle domain range of 45 ° to 135 °. The accuracy of the power estimation for a single signal source and 30 signal sources at different signal-to-noise ratios is shown in fig. 4. It can be seen that, under the two conditions of a single signal source and 30 signal sources, the power estimation error of the method provided by the invention is less than 1dB, and the estimation of the signal power can be accurately carried out.
Simulation example 2: the same co-prime array as the simulation example 1 is adopted to receive incident signals, the signal-to-noise ratio is set to be 10dB, the power estimation accuracy of a single signal source and 30 signal sources is tested under the condition of different sampling snapshot numbers, and other simulation conditions are kept the same as the simulation example 1. The power estimation accuracy for a single signal source and 30 signal sources under different sampling snapshot conditions is shown in fig. 5. It can be seen that the root mean square errors of the power estimates are all less than 1dB in the case of a single signal source; in the case of 30 signal sources, the root mean square error of the power estimate decreases below 1dB for fast sampling beats greater than 300. Therefore, the method provided by the invention can estimate the signal power more accurately.
Simulation example 3: the incident signal was received using the same relatively prime array as in simulation example 1, with the signal-to-noise ratio set to 10dB and the sample fast beat number set to T500. Assuming that the number of incident narrowband signals is 30, the incident directions are uniformly distributed in the spatial angle domain range of 45-135 deg. The spatial power spectrum obtained by the method of the present invention is shown in fig. 6, wherein the vertical dotted line represents the real direction of the incident signal source. It can be seen that the method provided by the present invention can effectively distinguish the 30 incident signal sources. For the traditional direction of arrival estimation method adopting a uniform linear array, 27 physical array elements are utilized to distinguish 26 incident signals at most, and the result shows that the method provided by the invention realizes the increase of the degree of freedom.
In summary, the method provided by the present invention obtains the spatial response by performing inverse discrete fourier transform on the second-order equivalent virtual signal described by the angle-spatial frequency, constructs a spatial power spectrum based on the spatial response, obtains the DOA estimation by performing a spectral peak search process on the constructed spatial power spectrum, and obtains the power estimation value of the corresponding signal while ensuring the DOA estimation freedom performance improvement. The operation of the inverse discrete Fourier transform effectively reduces the computational complexity, better meets the actual application scene with higher requirements on real-time performance, and is easy to realize in hardware in an actual system.
The above-described embodiments are intended to illustrate rather than to limit the invention, and any modifications and variations of the present invention are within the spirit of the invention and the scope of the appended claims.

Claims (7)

1. A direction of arrival estimation method based on inverse discrete Fourier transform of a second-order equivalent virtual signal of a co-prime array is characterized by comprising the following steps:
(1) the receiving end uses 2M + N-1 array elements and constructs according to a co-prime array structure; wherein M and N are relatively prime integers;
(2) suppose there are L from [ theta ═ theta1，θ2，…，θL]TDirectional far-field narrow-band incoherent signal source [. ]]TRepresenting a transposition operation, receiving an incident signal by using the co-prime array constructed in step (1), and modeling a co-prime array received signal x (t) at time t as follows:
Figure FDA0002533909660000011
wherein x (t) is a (2M + N-1) × 1-dimensional vector, sl(t) is the waveform of the first incident signal, n (t) is the noise component independent of each signal source, a (theta)l) To correspond to thetalThe guiding vector of the co-prime array of the directional signal source can be expressed as
Figure FDA0002533909660000012
Wherein, mui1,2,3, …,2M + N-1 represents the actual position of the ith physical array element in the co-prime array, and the position of the first physical array element is mu1λ is the wavelength of the incident narrowband signal, and j is the imaginary unit;
constructing a covariance matrix according to the co-prime array received signal x (t): sampling covariance matrix obtained using T sampling snapshots
Figure FDA0002533909660000013
For the theoretical covariance matrix RxCarrying out approximate substitution;
(3) deriving a second order equivalent virtual signal corresponding to the augmented virtual uniform linear array from the co-prime array received signal: sampling covariance matrix by vectorization
Figure FDA0002533909660000014
Obtaining a virtual array equivalent virtual signal y:
Figure FDA0002533909660000015
wherein,
Figure FDA0002533909660000016
is (2M + N-1)2The steering matrix is maintained at × L,
Figure FDA0002533909660000017
to include L incident signalsA vector of the power of the signal source,
Figure FDA0002533909660000018
for noise power, I ═ vec (I)2M+N-1)，I2M+N-1Representing a (2M + N-1) × (2M + N-1) dimensional identity matrix, vec (·) representing a vectorization operation, i.e., stacking columns in the matrix in sequence to form a new vector, (·)*It is meant a conjugate operation of the two,
Figure FDA0002533909660000021
represents the kronecker product;
non-uniform virtual array S corresponding to vector yDExpressed as:
SD＝{±(Mn-Nm)d，0≤n≤N-1，0≤m≤2M-1}，
where d is half the wavelength of the incident narrowband signal, i.e.
Figure FDA0002533909660000022
Selecting a non-uniform virtual array SDThe virtual array elements of the medium-maximum continuous part form a virtual uniform linear array S containing 2V +1 virtual array elementsV＝{-Vd，-(V-1)d，...，0，...，(V-1)d，Vd}，V＝MN+M-1；
Selecting the vector y corresponding to SVThe equivalent virtual signals of the positions of the virtual array elements form a second-order equivalent virtual signal z corresponding to the virtual uniform linear arrayθIt can be expressed as:
Figure FDA0002533909660000023
wherein B (θ) ═ B (θ)1)，b(θ2)，...b(θL)]Column l thereof
Figure FDA0002533909660000024
To correspond to thetalVirtual uniform linear array steering vector of directional signal source, eIs selected from i to correspond to SVVectors composed of elements of the medium array elements;
(4) defining an angle-space frequency and describing a second-order equivalent virtual signal of the virtual uniform linear array by the angle-space frequency; defining the angle-space frequency as the number of signal cycles of a narrow-band signal from a direction theta in space, which are propagated within the propagation distance difference range between adjacent array elements; second-order equivalent virtual signal z of the virtual uniform linear array in the step (3)θCan be equivalently expressed in the angular-spatial frequency domain as:
Figure FDA0002533909660000025
wherein, B (ξ) ═ B (ξ)1)，b(ξ2)，...b(ξL)]Column l thereof
Figure FDA0002533909660000026
Figure FDA0002533909660000027
ξ＝[ξ1，ξ2，…，ξL]TAngle-space frequencies corresponding to L angles contained in θ;
(5) performing inverse discrete Fourier transform on the second-order equivalent virtual signal described by the angle-space frequency, and constructing a space power spectrum: second order equivalent virtual signal z represented by angle-space frequency by inverse discrete Fourier transformξConverting the space domain to obtain a K × 1 dimensional space response psi;
constructing a spatial power spectrum, wherein the horizontal axis of the spectrum represents an angle theta, and the relation between the angle theta and the k-th element of the spatial response vector can be expressed as:
Figure FDA0002533909660000031
where K is 0, 1., K-1, K is the dimension of the spatial response ψ, arccos (·) is an inverse cosine function, and h is a guarantee
Figure FDA0002533909660000032
Coefficients satisfying the domain of the inverse cosine function when
Figure FDA0002533909660000033
When h is-1, when
Figure FDA0002533909660000034
When h is 0; the vertical axis of the spectrum represents the modulus p (k) of the kth element in the spatial response vector;
(6) estimating the direction of arrival and the signal power according to the obtained space power spectrum; and (5) performing spectrum peak searching operation on the spatial power spectrum in the step (5), wherein the angle corresponding to the first L peak values with the maximum amplitude is the direction of arrival estimation of the L incident signals, and the peak value amplitude is the power estimation value of the corresponding signal.
2. The method for estimating the direction of arrival of an inverse discrete fourier transform of a second-order equivalent virtual signal based on a co-prime array as claimed in claim 1, wherein: the coprime array structure in the step (1) can be specifically described as follows: firstly, selecting a group of coprime integers M, N to construct a pair of sparse uniform linear sub-arrays; the first subarray comprises 2M array elements with a spacing Nd, and the positions of the array elements are 0, Nd, (2M-1) Nd; the second subarray comprises N array elements with the interval Md, and the positions of the N array elements are 0, Md, (N-1) Md; and then combining the two sub-arrays according to a mode that the first array elements are overlapped to obtain a non-uniform co-prime array containing 2M + N-1 physical array elements.
3. The estimation method of the direction of arrival of the co-prime array based on the angle-space frequency domain fast fourier transform of claim 1, wherein: in step (3), if there are multiple different elements and S in yVIf the same virtual array element position corresponds to each other, any one of the elements is selected as a composition vector zθOf (2) is used.
4. The method for estimating the direction of arrival of an inverse discrete fourier transform of a second-order equivalent virtual signal based on a co-prime array as claimed in claim 1, wherein: the angle-space frequency in step (4) is a frequency related to the incident angle of the signal, and is specifically defined by the following way: the arrival of a narrow-band signal from the direction θ at an adjacent array element with a distance d in space generates a propagation distance difference u, which can be expressed as:
u＝dcosθ
under the condition that the array element spacing is fixed, the propagation distance difference u only changes along with the incident signal angle theta, so that the definition of the angle-space frequency in the step (4) is obtained; the angle-spatial frequency ξ versus incident signal angle θ can be expressed as:
Figure FDA0002533909660000041
5. the method for estimating the direction of arrival of an inverse discrete fourier transform of a second-order equivalent virtual signal based on a co-prime array as claimed in claim 1, wherein: the second-order equivalent virtual signal z obtained in the step (5) through inverse discrete Fourier transformξThe K × 1 dimensional spatial response ψ of can be expressed as:
Figure FDA0002533909660000042
wherein,
Figure FDA0002533909660000043
denotes the inverse discrete Fourier transform operation, FKIs a K-point inverse discrete fourier transform matrix, which can be expressed as:
Figure FDA0002533909660000044
the resulting spatial response ψ is a K × 1 dimensional vector, K being the dimension of the spatial response ψ.
6. The estimation method of the direction of arrival of the co-prime array based on the angle-space frequency domain fast fourier transform of claim 1, wherein: in the step (5), the constructed spatial spectrum reflects the response amplitude at each angle in space, where there are L peaks corresponding to L incident signals.
7. The method of estimating direction of arrival based on inverse discrete fourier transform of a second order equivalent virtual signal of a relatively prime array as claimed in claim 5, wherein: in the step (5), there are L peaks corresponding to L incident signals in the spatial power spectrum, and the conclusion is obtained as follows: the relation p (k) used to construct the spatial power spectrum may be specifically expressed as:
Figure FDA0002533909660000045
wherein,
Figure FDA0002533909660000046
for the power of the l signal, (-) represents an impulse function, and r is an integer and is used for representing that the impulse function is a periodic impulse string sequence; according to the nature of the impulse function, only when
Figure FDA0002533909660000047
When P (k) is a peak value; since the L incident signals have different angle-space frequencies, each of which can and can only cause one peak in the spatial power spectrum, there are L peaks in the spatial power spectrum corresponding to the L incident signals.
CN201810460228.7A 2018-05-15 2018-05-15 Direction-of-arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of co-prime array Active CN108710102B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810460228.7A CN108710102B (en) 2018-05-15 2018-05-15 Direction-of-arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of co-prime array

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810460228.7A CN108710102B (en) 2018-05-15 2018-05-15 Direction-of-arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of co-prime array

Publications (2)

Publication Number Publication Date
CN108710102A CN108710102A (en) 2018-10-26
CN108710102B true CN108710102B (en) 2020-09-04

Family

ID=63869008

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810460228.7A Active CN108710102B (en) 2018-05-15 2018-05-15 Direction-of-arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of co-prime array

Country Status (1)

Country Link
CN (1) CN108710102B (en)

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109471065B (en) * 2018-09-28 2021-05-25 中国电子科技集团公司第三十六研究所 Direction finding method for coherent signals
CN110133574B (en) * 2019-07-02 2022-12-16 华南理工大学 One-dimensional DOA estimation method utilizing secondary virtual expansion of multi-frequency signals
CN111610486B (en) * 2020-05-03 2022-12-27 浙江大学 High-resolution accurate two-dimensional direction of arrival estimation method based on planar co-prime array virtual domain tensor spatial spectrum search
WO2021068495A1 (en) 2020-05-03 2021-04-15 浙江大学 Degree of freedom-enhanced spatial spectrum estimation method based on block sampled tensor signal construction using planar co-prime array
CN112162240A (en) * 2020-09-25 2021-01-01 五邑大学 Sparse frequency waveform generation method and device based on co-prime array and storage medium
CN114167354B (en) * 2021-11-25 2025-01-14 厦门大学 A method for underwater DOA estimation based on high degree of freedom
CN114444298A (en) * 2022-01-21 2022-05-06 浙江大学 Estimation method for two-dimensional direction of arrival of co-prime area array based on virtual domain tensor filling
CN115407128B (en) * 2022-07-29 2023-10-27 南京航空航天大学 Electric power system harmonic wave and inter-harmonic wave frequency estimation method based on inter-mass sampling

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107015190A (en) * 2017-03-01 2017-08-04 浙江大学 Relatively prime array Wave arrival direction estimating method based on the sparse reconstruction of virtual array covariance matrix
CN107422295A (en) * 2017-08-30 2017-12-01 浙江大学 The Wave arrival direction estimating method represented based on relatively prime array virtual Domain equivalent signal atom norm
CN108680892A (en) * 2018-05-15 2018-10-19 浙江大学 Relatively prime array Wave arrival direction estimating method based on angle-spatial frequency domain Fast Fourier Transform (FFT)

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107015190A (en) * 2017-03-01 2017-08-04 浙江大学 Relatively prime array Wave arrival direction estimating method based on the sparse reconstruction of virtual array covariance matrix
CN107422295A (en) * 2017-08-30 2017-12-01 浙江大学 The Wave arrival direction estimating method represented based on relatively prime array virtual Domain equivalent signal atom norm
CN108680892A (en) * 2018-05-15 2018-10-19 浙江大学 Relatively prime array Wave arrival direction estimating method based on angle-spatial frequency domain Fast Fourier Transform (FFT)

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
Doa estimation by covariance matrix sparse reconstruction of coprime array;Chengwei Zhou等;《 2015 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)》;20150806;2369-2373 *
High-resolution Direction-of-Arrival estimation in SNR and snapshot challenged scenarios using multi-frequency coprime arrays;Yang Liu等;《 2017 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)》;20170619;3434-3438 *
Sparsity-based DOA estimation using co-prime arrays;Yimin D. Zhang等;《 2013 IEEE International Conference on Acoustics, Speech and Signal Processing》;20131221;3967-3971 *

Also Published As

Publication number Publication date
CN108710102A (en) 2018-10-26

Similar Documents

Publication Publication Date Title
CN108710102B (en) Direction-of-arrival estimation method based on second-order equivalent virtual signal inverse discrete Fourier transform of co-prime array
CN107315160B (en) Relatively prime array Wave arrival direction estimating method based on interpolation virtual array signal atom norm minimum
CN107329108B (en) The relatively prime array Wave arrival direction estimating method rebuild based on interpolation virtual array covariance matrix Toeplitzization
CN106324558B (en) Broadband signal DOA estimation method based on co-prime array
CN107015190A (en) Relatively prime array Wave arrival direction estimating method based on the sparse reconstruction of virtual array covariance matrix
CN107290709B (en) The relatively prime array Wave arrival direction estimating method decomposed based on vandermonde
CN107561484B (en) Direction-of-arrival estimation method based on interpolation co-prime array covariance matrix reconstruction
CN107092004A (en) Relatively prime array Wave arrival direction estimating method based on signal subspace rotational invariance
CN113189592B (en) Vehicle-mounted millimeter wave MIMO radar angle measurement method considering amplitude mutual coupling error
CN107450047A (en) Compressed sensing DOA estimation method based on unknown mutual coupling information under nested battle array
CN107037392A (en) A kind of relatively prime array Wave arrival direction estimating method of free degree increase type based on compressed sensing
CN110133574B (en) One-dimensional DOA estimation method utilizing secondary virtual expansion of multi-frequency signals
CN113835063B (en) Unmanned aerial vehicle array amplitude and phase error and signal DOA joint estimation method
CN107104720A (en) The relatively prime array adaptive beamforming method rebuild based on covariance matrix virtual Domain discretization
CN113296049A (en) Conjugate amplification DOA estimation method of non-circular signals in co-prime array pulse environment
CN111366893A (en) Non-circular signal azimuth angle estimation method under unknown mutual coupling condition of uniform circular array
CN116299150B (en) Two-dimensional DOA estimation method of dimension-reduction propagation operator in uniform area array
CN111352063A (en) Two-dimensional direction finding estimation method based on polynomial root finding in uniform area array
CN108614234B (en) Direction-of-arrival estimation method based on multi-sampling snapshot co-prime array received signal fast Fourier inverse transformation
CN115236589B (en) Polar region under-ice DOA estimation method based on covariance matrix correction
CN113791379B (en) Orthogonal matching pursuit DOA estimation method under nested array non-Gaussian environment
CN109471087B (en) Direction-of-arrival estimation method based on co-prime MIMO radar difference set and signal collection fast Fourier transform
CN108680892B (en) Estimation method of direction of arrival of co-prime array based on angle-space frequency domain fast Fourier transform
CN109298381A (en) A kind of relatively prime battle array coherent signal azimuth estimation method based on variational Bayesian
CN114563760B (en) Second-order super-beam forming method, equipment and medium based on SCA array

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant