CN114755628A - Method for estimating direction of arrival of acoustic vector sensor array under non-uniform noise - Google Patents

Method for estimating direction of arrival of acoustic vector sensor array under non-uniform noise Download PDF

Info

Publication number
CN114755628A
CN114755628A CN202210353033.9A CN202210353033A CN114755628A CN 114755628 A CN114755628 A CN 114755628A CN 202210353033 A CN202210353033 A CN 202210353033A CN 114755628 A CN114755628 A CN 114755628A
Authority
CN
China
Prior art keywords
covariance matrix
signal
sparse
noise
vector
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.)
Pending
Application number
CN202210353033.9A
Other languages
Chinese (zh)
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.)
Henan University of Technology
Original Assignee
Henan University of Technology
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 Henan University of Technology filed Critical Henan University of Technology
Priority to CN202210353033.9A priority Critical patent/CN114755628A/en
Publication of CN114755628A publication Critical patent/CN114755628A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Direction-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/80Direction-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/8003Diversity systems specially adapted for direction finding
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/16Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/18Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Computational Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Software Systems (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Computing Systems (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Operations Research (AREA)
  • Probability & Statistics with Applications (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Abstract

The invention relates to the technical field of array signal processing, and discloses a method for estimating the direction of arrival of an acoustic vector sensor array under non-uniform noise, which comprises the following steps: step 1, receiving an output vector X output by an acoustic vector sensor array; step 2, based on a sparse signal model, calculating a signal covariance matrix P and a noise covariance matrix Q according to the output vector X, and constructing a sparse covariance matrix R by using the obtained signal covariance matrix P and the noise covariance matrix Q; step 3, constructing a cost function based on the sparse covariance matrix R, estimating sparse signal power by using the cost function, and obtaining a signal power vector
Figure DDA0003581551110000012
Step 4, for signal power vector
Figure DDA0003581551110000011
And searching a spectral peak, wherein the sound source position corresponding to the obtained spectral peak is the estimated target azimuth. The method can improve the performance of the orientation estimation of the acoustic vector sensor array under the nonuniform noise.

Description

Method for estimating direction of arrival of acoustic vector sensor array under non-uniform noise
Technical Field
The invention relates to the technical field of array signal processing, in particular to a method for estimating the direction of arrival of an acoustic vector sensor array under non-uniform noise.
Background
The acoustic vector sensor is a novel sensor, and can be used for measuring sound pressure and vibration velocity information in a sound field synchronously at the same point in space. Compared with the traditional sound pressure sensor, the sound vector sensor can acquire sound pressure information in a sound field and sound vibration velocity information, and has better signal detection and direction estimation performance by utilizing the additional sound vibration velocity information. Compared with a single acoustic vector sensor, the acoustic vector sensor array formed by the acoustic vector sensor has the advantages of high signal gain, strong interference suppression capability, high spatial resolution and the like. The method is widely applied to the detection fields of sonar, radar and the like.
Estimating the Direction Of Arrival (DOA) Of an underwater sound source is a popular research Direction in the field Of passive detection. The beam forming method is a more classical DOA estimation method, has the advantages of simple structure and small calculation amount, and is limited by physical aperture and has poor resolution performance. Although a subspace-based algorithm, such as a Multiple Signal Classification (MUSIC) algorithm and its improved algorithm, can achieve high-resolution estimation, the DOA estimation performance of the method is further deteriorated in low snr, few snapshots, and coherent source scenarios. The sparse reconstruction method utilizes the characteristic that incident signals have sparse distribution relative to the whole space to perform grid division on the angle space where the signals are located, and converts the azimuth estimation problem into the sparse signal reconstruction problem. The algorithm makes up the defect that the azimuth estimation performance is deteriorated under the conditions of coherent source, low signal-to-noise ratio and less snapshots in the traditional azimuth estimation method. l 1The norm optimization method is a relatively classical sparse signal reconstruction method and has higher estimation than the traditional method under the conditions of low signal-to-noise ratio and less snapshotsMeter accuracy and robustness, however l1The norm optimization method faces the difficult problem that regularization parameters are difficult to select properly in a complex underwater acoustic environment. An Iterative Adaptive Approach (IAA) is a sparse method without setting user parameters, has the advantages of a sparse reconstruction algorithm and a smaller calculation amount, and is widely applied to radar systems. However, the above method is based on an ideal assumption that the noise powers of the sound pressure and the sound vibration velocity components are equal to each other, and the azimuth estimation performance of the acoustic vector sensor array is studied. Because the noise received by the acoustic vector sensor is mainly environmental noise rather than noise inside the sensor, and the vibration speed axis of the acoustic vector sensor has directivity, part of the environmental noise measured in the direction of the vibration speed axis is filtered, so that the noise power of the vibration speed component is not equal to the noise power of the sound pressure component. In this case, the azimuth estimating performance of the above method will be seriously deteriorated.
To solve the problem of degraded performance of sound pressure sensor array orientation Estimation under non-uniform Noise, Bin Liao et al propose an Iterative Subspace Estimation method in the paper "Iterative Methods for Subspace and DOA Estimation in Noise" (IEEE Transactions on Signal Processing Volume 64, Issue 12.2016.PP 3008-. Although the method can be applied to the acoustic vector sensor array to solve the problem of deteriorated azimuth estimation performance of the acoustic vector sensor array under the condition of non-uniform noise, the two methods adopt an iterative mode to estimate the noise covariance matrix, so that the calculation complexity of the algorithm is high, and the requirement of the system on the real-time property is difficult to meet. In order to solve the problem, in the article "Augmented subspace MUSIC method for DOA estimation using an acoustic vector sensor array" (IET Radar, sound & Navigation Volume 13, Issue 6.2019.PP 969-. However, this algorithm is only applicable to the case where the noise power is different between the sound pressure and the vibration velocity in the large aperture array and the acoustic vector sensor. When the sensors are correspondingly different, the sound pressure and the vibration velocity between the acoustic vector sensors and the noise received between the vibration velocity and the vibration velocity channels are different, and in this case, the azimuth estimation performance of the asmuscic method is also seriously deteriorated.
Disclosure of Invention
In view of the above situation, an object of the present invention is to provide a regularized Sparse Covariance Matrix Fitting (WSCMF) orientation estimation method, so as to improve the accuracy of orientation estimation of an acoustic vector sensor array under non-uniform noise.
In order to achieve the purpose, the invention is realized by the following technical scheme:
a method for estimating the direction of arrival of an acoustic vector sensor array under non-uniform noise comprises the following steps:
step 1, receiving an output vector X output by an acoustic vector sensor array;
step 2, based on a sparse signal model, calculating a signal covariance matrix P and a noise covariance matrix Q according to the output vector X, and constructing a sparse covariance matrix R by using the obtained signal covariance matrix P and the noise covariance matrix Q;
step 3, constructing a cost function based on the sparse covariance matrix R, estimating sparse signal power by using the cost function, and obtaining a signal power vector
Figure BDA0003581551090000041
Step 4, for signal power vector
Figure BDA0003581551090000042
And searching a spectral peak, wherein the sound source position corresponding to the obtained spectral peak is the estimated target azimuth.
Further, in the above method for estimating the direction of arrival of the acoustic vector sensor array under the non-uniform noise, step 2 includes the following substeps:
Step 21, constructing a sparse signal model of the acoustic vector sensor array under the condition of non-uniform noise:
supposing that K far-field equipower narrowband signals with the wavelength of lambda are incident on a uniform linear array consisting of M acoustic vector sensors, the spacing between array elements is d, and the vibration velocity components of the narrowband signals along the directions of an x axis and a y axis are u respectivelyxAnd uyAnd the X-axis and the y-axis have orthogonality, the output vector X of the acoustic vector sensor array can be expressed as L in the case of fast beat number
X=A(θ)S+W
Wherein X is [ X (1), X (2), …, X (L) ]]Represents an output vector of the acoustic vector sensor array, and a (θ) ═ a (θ)1),a(θ2),…,a(θK)]An array manifold matrix representing an array of ideal acoustic vector sensors, θ ═ θ [ θ ]1,θ2,…,θK]The azimuth angle of the target is represented,
Figure BDA0003581551090000043
Figure BDA0003581551090000044
a manifold vector representing the kth target,
Figure BDA0003581551090000045
representing the Cartesian product, h (θ)k)=[1 cosθk sinθk]TAn array manifold vector representing the incidence of the kth target on the acoustic vector transducer,
Figure BDA0003581551090000046
Figure BDA0003581551090000047
denotes a sound pressure response coefficient of the kth target incident on the mth acoustic vector sensor, W ═ W (1), W (2), …, W (l)]Represents a noise vector, S ═ S (1), S (2), …, S (l)]Representing a vector of input signal waveforms;
it is assumed that the angular space (-180 deg.) is uniformly divided into NA discrete grid is formed, and the number N of the discrete grid is far larger than the number K of the sound sources, the incoming wave directions of all the sound sources can be expressed as
Figure BDA0003581551090000051
The sparse signal model of the acoustic vector sensor array may be represented as
Figure BDA0003581551090000052
In the formula,
Figure BDA0003581551090000053
an array manifold matrix representing an array of acoustic vector sensors under a sparse signal model,
Figure BDA0003581551090000054
is a row sparse matrix and its non-zero rows represent the locations of the real signals;
step 22, calculating a signal covariance matrix P and a noise covariance matrix Q:
the signal covariance matrix P can be expressed as
Figure BDA0003581551090000055
Wherein E {. denotes the desired operation, diag {. denotes the diagonal matrix operation,
Figure BDA0003581551090000056
representing the power of the corresponding signal of the nth grid;
the noise covariance matrix Q can be expressed as
Q=diag{q1,…,qM}
In the formula,
Figure BDA0003581551090000057
assuming that the incident signals are uncorrelated with each other, the covariance matrix of the acoustic vector sensor array under the sparse signal model can be generally expressed as
Figure BDA0003581551090000058
In general, if the background noise is white gaussian noise, the noise covariance matrix Q is σ2I3M,σ2Representing the noise power, I, output by each channel of the acoustic vector sensor array in a uniform noise environment3MRepresenting a 3M x 3M identity matrix. However, in practical applications, since the noise received by the acoustic vector sensor is mainly environmental noise rather than noise inside the sensor, and the vibration speed axis of the acoustic vector sensor has directivity, a part of the environmental noise measured in the vibration speed axis direction is filtered, so that the noise power of the acoustic vibration speed component is not equal to the noise power of the sound pressure component. In addition, the received noise power between the sound pressure channel and the sound velocity channel between the acoustic vector sensors is made different in consideration of the difference in response between the sensors, that is, the noise power is different between the sound pressure channel and the sound velocity channel
Figure BDA0003581551090000061
In order to solve the problem, the invention defines a new array manifold matrix of the acoustic vector sensor array, and reconstructs a sparse covariance matrix R, so that the sparse covariance matrix R can simultaneously estimate the signal power on each grid and the noise power output by each channel in the acoustic vector sensor array. It should be noted that each channel of the acoustic vector sensor array refers to a sound pressure channel of each acoustic vector sensor array element, and a sound vibration velocity channel corresponding to a sound vibration velocity component of each acoustic vector sensor array element along the x-axis and the y-axis directions, respectively.
Step 23, constructing a sparse covariance matrix R:
in order to estimate the power of sparse signals and the noise power output by each channel in an acoustic vector sensor array, a new array manifold matrix of the acoustic vector sensor array is defined
Figure BDA0003581551090000062
In which I3MRepresenting a 3M x 3M identity matrix, assuming that the incident signals are uncorrelated with each other, the sparse covariance matrix R can be represented as
R=DΞDH
Wherein,
Figure BDA0003581551090000063
Figure BDA0003581551090000071
in the formula, the first N term of the XI diagonal element represents the signal power, and the N +1 to N +3M terms represent the noise power output by each channel in the acoustic vector sensor array.
Thus, the noise covariance matrix Q can be further expressed as
Figure BDA0003581551090000072
In the formula DnColumn n data representing D, DN+l=ul,l=1,...,3M,ulIs represented by I3MColumn l.
Further, in the above method for estimating the direction of arrival of the acoustic vector sensor array under the non-uniform noise, step 3 includes the following substeps:
step 31, constructing an interference-plus-noise covariance matrix B according to the sparse covariance matrix Rn
According to the sparse covariance matrix R constructed in step 23, the signal covariance matrix corresponding to the nth grid can be expressed as
Figure BDA0003581551090000073
The interference plus noise covariance matrix corresponding to the signal on the nth grid is Bn=R-Rn=DΞDH-Rn
Step 32, constructing the following cost function based on the fitting criterion of the sparse covariance matrix:
Figure BDA0003581551090000081
in the formula (I)
Figure BDA0003581551090000082
Representing the covariance matrix fitting term, the second term
Figure BDA0003581551090000083
Represents a sparse signal compensation term, | · | non-wovenFIt means that the F-norm operation is taken,
Figure BDA0003581551090000084
a sample covariance matrix is represented as a function of,
Figure BDA0003581551090000085
q denotes a user parameter, λ, for compensating the signalsRegularization parameter representing the relationship of sparsity of the equilibrium signal to fitting error, ands>0;
step 33, solving the sparse signal power by using an iterative method, and obtaining a signal power vector
Figure BDA0003581551090000086
Further, in the above method for estimating the direction of arrival of the acoustic vector sensor array under non-uniform noise, step 33 includes the following substeps:
step 331, in order to reduce the influence of the computation complexity and noise on the covariance matrix fitting term effect in the cost function, the covariance matrix of the sample is subjected to
Figure BDA0003581551090000087
Singular value decomposition is carried out to obtain
Figure BDA0003581551090000088
In the formula,
Figure BDA0003581551090000089
a signal sub-space is represented that is,
Figure BDA00035815510900000810
representing a noise subspace;
the cost function can be expressed as
Figure BDA00035815510900000811
Step 332, solving the sparse signal power by using an iterative method, setting the j +1 th iteration to obtain a unique solution, and obtaining the sparse signal power
Figure BDA00035815510900000812
In the formula, anIs composed of
Figure BDA00035815510900000813
In the form of a short-hand writing of (1),
Figure BDA00035815510900000814
Figure BDA00035815510900000815
Figure BDA00035815510900000816
Figure BDA0003581551090000091
specifically, the process of solving the sparse signal power by using the iterative method is as follows:
as can be seen from the expression of the cost function, the introduction of the user parameter q makes the cost function about the variable to be solved
Figure BDA0003581551090000092
Is a non-linear function of (a). Therefore, solving directly from the cost function
Figure BDA0003581551090000093
Is very difficult. To solve this problem, assume that f (x) xqIs a concave function with respect to variable x. When the variable x > 0, for x0Is greater than 0, can be obtained
f(x)≤f(x0)+f′(x0)(x-x0) (1)
Order to
Figure BDA0003581551090000094
The superscript j represents the number of iterations, and can be substituted for the formula (1)
Figure BDA0003581551090000095
Multiplying both sides of equation (2) by λsAnd then the covariance matrix fitting term
Figure BDA0003581551090000096
Can be further converted into
Figure BDA0003581551090000097
In the formula,
Figure BDA0003581551090000098
thus, the cost function is related to
Figure BDA0003581551090000099
The minimized problem is ultimately equivalent to
Figure BDA00035815510900000910
In the formula (4), the reaction mixture is,
Figure BDA00035815510900000911
is about
Figure BDA00035815510900000912
A function of
Figure BDA00035815510900000913
In a clear view of the above, it is known that,
Figure BDA00035815510900000914
also relates to
Figure BDA00035815510900000915
As a function of (c). Therefore, it can be considered that in the formula (4),
Figure BDA00035815510900000916
is an unknown variable. It is assumed that in the formula (4),
Figure BDA0003581551090000101
is estimated in the j-th iteration, equation (5) is minimized
Figure BDA0003581551090000102
Can be converted into minimization
Figure BDA0003581551090000103
In the formula,
Figure BDA0003581551090000104
as is apparent from the formula (7),
Figure BDA0003581551090000105
is about
Figure BDA0003581551090000106
Is a linear function of (a). Therefore, can be directly selected from
Figure BDA0003581551090000107
The sparse signal power and the noise power are solved.
The formula (7) is further simplified according to the property of Frobenius norm
Figure BDA0003581551090000108
In the formula, Tr {. cndot } represents a matrix tracing operation. To pair
Figure BDA0003581551090000109
About
Figure BDA00035815510900001010
The first derivative is obtained
Figure BDA00035815510900001011
To pair
Figure BDA00035815510900001012
About
Figure BDA00035815510900001013
The second derivative is obtained
Figure BDA00035815510900001014
It can be obtained from the formula (10),
Figure BDA00035815510900001015
this is always true. Therefore, the temperature of the molten metal is controlled,
Figure BDA00035815510900001016
there is a unique solution. Order to
Figure BDA0003581551090000111
Available function
Figure BDA0003581551090000112
The only solution in the j +1 th iteration is
Figure BDA0003581551090000113
As can be seen from equation (11), when calculating the corresponding sparse signal power on N discrete grids, it is necessary to calculate the interference-plus-noise covariance matrix B N timesnThe inverse of (2) greatly increases the amount of computation. To reduce the amount of calculation, Woodbury's formula can be used
Figure BDA0003581551090000114
By substituting formula (12) for formula (11), a compound of formula (I) can be obtained
Figure BDA0003581551090000115
In the formula,
Figure BDA0003581551090000116
Figure BDA0003581551090000117
Figure BDA0003581551090000118
the sparse signal power can be obtained according to the definition of D
Figure BDA0003581551090000119
In the formula, anIs composed of
Figure BDA00035815510900001110
In a simplified form.
Then the noise power
Figure BDA00035815510900001111
In the formula,
Figure BDA0003581551090000121
Figure BDA0003581551090000122
thus, the sound pressure channel and the sound velocity channel have a noise power of
Figure BDA0003581551090000123
Figure BDA0003581551090000124
Figure BDA0003581551090000125
Step 333, obtaining the signal power vector of the signals corresponding to the N grids according to the sparse signal power obtained in step 332 as
Figure BDA0003581551090000126
Compared with the prior art, the invention has the beneficial effects that:
1. The method defines a new array manifold matrix of the acoustic vector sensor array, reconstructs a sparse covariance matrix containing sparse signal power and noise power according to the newly defined array manifold matrix of the acoustic vector sensor array, constructs a cost function about the sparse signal power and the noise power based on a sparse covariance matrix fitting criterion and a sparse signal compensation mechanism, and estimates the noise power and the signal power by adopting a regularized weighted covariance matrix fitting method. In each iteration, the sparse covariance matrix is reconstructed through the sparse signal power and the noise power estimated by the method provided by the invention, so that the sparse signal power and the noise power estimated in the next iteration are more accurate, and the sparse signal power is subjected to spectral peak search when the iteration is terminated, so that more accurate azimuth estimation can be realized.
Compared with the method for estimating the signal power only in the existing direction of arrival estimation method, the method for estimating the signal power based on the noise power has the advantages that the signal power estimated by utilizing the accurate noise power is more accurate than that of the existing estimation method by reconstructing the sparse covariance matrix containing the sparse signal power and the noise power and estimating the signal power and the noise power at the same time. .
2. Compared with the whitening processing of the output vector of the acoustic vector sensor array in the existing direction of arrival estimation method, the direction of arrival estimation method provided by the invention has less calculation amount and higher estimation precision. In practical application, the method provided by the invention can improve the orientation estimation precision of the acoustic vector sensor array under the nonuniform noise.
Drawings
FIG. 1 is a flow chart of a WSCMF method according to the present invention;
FIG. 2 is a uniform linear array model of an acoustic vector sensor array;
FIG. 3 is a normalized space spectrogram of the WSCMF method and the MUSIC method, the IAA method, the IMLSE method and the ASMUSIC method proposed by the present invention;
FIG. 4 is a graph showing the relationship between the root mean square error and the angular interval of the WSCMF method, the MUSIC method, the IAA method, the IMLSE method and the ASMUSIC method according to the present invention;
FIG. 5 is a graph showing the relationship between the root mean square error and the signal-to-noise ratio of the WSCMF method, the MUSIC method, the IAA method, the IMLSE method and the ASMUSIC method according to the present invention;
FIG. 6 is a graph showing the relationship between the root mean square error and the worst noise power ratio of the WSCMF method, the MUSIC method, the IAA method, the IMLSE method and the ASMUSIC method according to the present invention.
Detailed Description
The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be obtained by a person skilled in the art without making any creative effort based on the embodiments in the present invention, belong to the protection scope of the present invention.
Examples
As shown in fig. 1-6, a method for estimating the direction of arrival of an acoustic vector sensor array under non-uniform noise comprises the following steps:
step 1, receiving an output vector X output by an acoustic vector sensor array;
step 2, based on a sparse signal model, calculating a signal covariance matrix P and a noise covariance matrix Q according to the output vector X, and constructing a sparse covariance matrix R by using the obtained signal covariance matrix P and the noise covariance matrix Q;
specifically, step 2 includes the following substeps:
step 21, constructing a sparse signal model of the acoustic vector sensor array under the condition of non-uniform noise:
assuming that K far-field equal-power narrow-band signals with the wavelength of lambda are incident on a uniform linear array consisting of M acoustic vector sensors, the arrangement of the acoustic vector sensor array is shown in FIG. 2, wherein the array element spacing is d, and the narrow-band signals are along The vibration velocity components in the x-axis and y-axis directions are respectively uxAnd uyAnd the X-axis and the y-axis have orthogonality, the output vector X of the acoustic vector sensor array can be expressed as L in the case of fast beat number
X=A(θ)S+W
Wherein X is [ X (1), X (2), …, X (L)]Representing an output vector of the acoustic vector sensor array, a (θ) ═ a (θ)1),a(θ2),…,a(θK)]Array manifold matrix representing an ideal acoustic vector sensor array, θ ═ θ1,θ2,…,θK]The azimuth of the target is represented and,
Figure BDA0003581551090000141
Figure BDA0003581551090000142
a manifold vector representing the kth target,
Figure BDA0003581551090000143
representing the Cartesian product, h (θ)k)=[1 cosθk sinθk]TAn array manifold vector representing the incidence of the kth target on the acoustic vector transducer,
Figure BDA0003581551090000151
Figure BDA0003581551090000152
denotes a sound pressure response coefficient of the kth target incident on the mth acoustic vector sensor, W ═ W (1), W (2), …, W (l)]Represents a noise vector, S ═ S (1), S (2), …, S (l)]Representing a vector of input signal waveforms;
assuming that the angle space (-180 deg.) is uniformly divided into N discrete grids, and the number N of the discrete grids is far larger than the number K of the sound sources, the incoming wave directions of all the sound sources can be expressed as
Figure BDA0003581551090000153
The sparse signal model of the acoustic vector sensor array may be represented as
Figure BDA0003581551090000154
In the formula,
Figure BDA0003581551090000155
an array manifold matrix representing an array of acoustic vector sensors under a sparse signal model,
Figure BDA0003581551090000156
is a row sparse matrix and its non-zero rows represent the locations of the real signals;
Step 22, calculating a signal covariance matrix P and a noise covariance matrix Q:
the signal covariance matrix P can be expressed as
Figure BDA0003581551090000157
Wherein E {. denotes the desired operation, diag {. denotes the diagonal matrix operation,
Figure BDA0003581551090000158
representing the power of the corresponding signal of the nth grid;
the noise covariance matrix Q can be expressed as
Q=diag{q1,…,qM}
In the formula,
Figure BDA0003581551090000159
step 23, constructing a sparse covariance matrix R:
array manifold matrix defining a new acoustic vector transducer array
Figure BDA00035815510900001510
Wherein I3MRepresenting 3M by 3MIdentity matrix, then sparse covariance matrix R can be expressed as
R=DΞDH
Wherein,
Figure BDA0003581551090000161
in the formula, the first N terms of the diagonal elements of xi represent the signal power, and the N +1 to N +3M represent the noise power output by each channel in the acoustic vector sensor array. Thus, the noise covariance matrix Q can be further expressed as
Figure BDA0003581551090000162
In the formula, DnColumn n data, DN, representing D+l=ul,l=1,...,3M,ulIs represented by3MColumn l.
Step 3, constructing a cost function based on the sparse covariance matrix R, estimating sparse signal power by using the cost function, and obtaining a signal power vector
Figure BDA0003581551090000163
Specifically, step 3 includes the following substeps:
step 31, an interference plus noise covariance matrix B is constructed according to the sparse covariance matrix Rn
According to the sparse covariance matrix R constructed in step 23, the signal covariance matrix corresponding to the nth grid can be expressed as
Figure BDA0003581551090000171
The interference plus noise covariance matrix corresponding to the signal on the nth grid is Bn=R-Rn=DΞDH-Rn
Step 32, based on the fitting criterion of the sparse covariance matrix, constructing the following cost function:
Figure BDA0003581551090000172
in the formula, | \ | non-countingFThe operation of taking the F norm is shown,
Figure BDA0003581551090000173
representing the sample covariance matrix, q representing the user parameter for the compensation signal, lambdasRegularization parameter representing the relationship of sparsity of the equilibrium signal to fitting error, and λs>0;
Step 33, solving the sparse signal power by using an iterative method, and obtaining a signal power vector
Figure BDA0003581551090000174
Specifically, step 33 includes the following substeps:
step 331, singular value decomposition is carried out on the sample covariance matrix to obtain
Figure BDA0003581551090000175
In the formula,
Figure BDA0003581551090000176
a signal sub-space is represented that is,
Figure BDA0003581551090000177
representing a noise subspace;
the cost function can be expressed as
Figure BDA0003581551090000178
Step 332, solving the sparse signal power by using an iterative method, setting the j +1 th iteration to obtain a unique solution, and obtaining the sparse signal power
Figure BDA0003581551090000179
In the formula, anIs composed of
Figure BDA00035815510900001710
In the form of a short-hand writing of (1),
Figure BDA0003581551090000181
Figure BDA0003581551090000182
Figure BDA0003581551090000183
step 333, obtaining the signal power vector of the signals corresponding to the N grids according to the sparse signal power obtained in step 332 as
Figure BDA0003581551090000184
Step 4, for signal power vector
Figure BDA0003581551090000185
And searching a spectral peak, wherein the sound source position corresponding to the obtained spectral peak is the estimated target azimuth.
Compared with the existing non-uniform noise power estimation method, the method designed by the invention has less calculation amount and higher estimation precision. Meanwhile, the signal power estimated by using the accurate noise power is more accurate than that of the existing estimation method.
The effect of the invention can be compared and explained with the existing MUSIC method, IAA method, Iterative Maximum Likelihood Subspace Estimation (IMLSE) method and ASMUSIC method through simulation experiment:
comparative example 1
The resolution performance of the WSCMF method provided by the invention is compared with that of the existing MUSIC method, IAA method, IMLSE method and ASMUSIC method on the normalized space spectrum.
In the present comparative example, the acoustic vector sensor array was set to a uniform linear array consisting of 4 acoustic vector sensors, the background noise of the acoustic vector sensor array was non-uniform noise, and the noise covariance matrix was set to
Q=diag{3,1.5,2,5,2,1.5,4,2,5,2,30,2,1.5}
Therefore, the Worst Noise Power Ratio (WNPR) is
Figure BDA0003581551090000191
In the formula,
Figure BDA0003581551090000192
the maximum value of the Q is represented by,
Figure BDA0003581551090000193
represents the minimum value in Q. The signal-to-noise ratio (SNR) is defined as
Figure BDA0003581551090000194
In the formula,
Figure BDA0003581551090000195
representing the signal power.
In the present comparative example, three far-field narrow-band signals with a signal-to-noise ratio of 10dB were set from the azimuth angle θ, respectively1=-15°,θ2=38°,θ3The sound velocity is 1500m/s, the system sampling frequency is 5kHz, the central frequency of the signal is 700Hz, and the fast beat number is 300. In the WSCMF method provided by the invention, the iteration number is 30, and the convergence threshold is 10 -3
The resolution performance of the WSCMF method provided by the invention and the existing MUSIC method, IAA method, IMLSE method and ASMUSIC method are compared on the normalized spatial spectrum, and the simulation experiment result is shown in figure 3.
From the normalized spatial spectrum shown in FIG. 3, it can be seen that neither the MUSIC method, the IAA method, or the IMLSE method can distinguish the location at θ238 deg. of target orientation. Although the ASMUSIC method can recognize that the azimuth is at θ2The target is at 38 °, but the estimated result deviates significantly from the true target azimuth angle. The WSCMF method provided by the invention can distinguish three targets, and the estimation result is closest to the azimuth angle of the real target, thereby showing that the method has better distinguishing performance.
Comparative example 2
The Root Mean Square Error (RMSE) is used as a criterion for judging the azimuth estimation performance, and the estimation accuracy of the WSCMF method provided by the invention is compared with the estimation accuracy of the existing MUSIC method, IAA method, IMLSE method and ASMUSIC method.
In general, a smaller root mean square error indicates a more accurate estimate, and a larger root mean square error indicates a less accurate estimate. The root mean square error can be expressed as
Figure BDA0003581551090000201
In the formula, thetakFor the true horizontal azimuth angle of the kth target,
Figure BDA0003581551090000202
The estimated horizontal azimuth angle of the kth target in the U-th monte carlo experiment, U, is the number of independent monte carlo experiments performed in the 3000 experiment.
In this comparative example, the estimation performance of several azimuth estimation methods was studied from the relationship curve of root mean square error and angle interval, the relationship curve of root mean square error and signal-to-noise ratio, and the relationship curve of root mean square error and worst noise power ratio, respectively, and the simulation experimental conditions were specifically set as follows:
(1) two far-field narrow bands with equal power and 5dB signal-to-noise ratioThe signals are incident on the acoustic vector sensor array, and the angles of the two incident signals are theta1=-14°,θ2The relationship between the root mean square error and the angle interval of the estimation results obtained by applying the WSCMF method provided by the invention and the existing MUSIC method, IAA method, IMLSE method and asmusch method is compared, and the result is shown in fig. 4;
(2) two far-field narrow-band signals with equal power are respectively from theta1=-14°,θ2When the signal-to-noise ratio is changed from-15 dB to 15dB when the signal-to-noise ratio is incident on the acoustic vector sensor array at 38 degrees, the relation between the root mean square error and the signal-to-noise ratio of the estimation result obtained by applying the WSCMF method provided by the invention and the existing MUSIC method, IAA method, IMLSE method and ASMUSIC method is compared, and the result is shown in FIG. 5;
(3) Two far-field narrow-band signals with equal power are respectively from the azimuth angle theta1=-14°,θ2When the signal-to-noise ratio is 5dB and the worst noise power ratio is changed from 10 to 80 in increments of 10 when the signal is incident on the acoustic vector sensor array at 38 °, the relationship between the root mean square error and the worst noise power ratio of the estimation results obtained by applying the WSCMF method provided by the present invention and the existing MUSIC method, IAA method, IMLSE method and asmusch method is compared, and the result is shown in fig. 6.
As can be seen from fig. 4, the IMLSE method and the asmusch method cannot distinguish two targets whose azimuth angle intervals are less than 45 ° and 35 °, respectively. The IAA method cannot resolve two targets with an angular separation of less than 30 °. The MUSIC method cannot resolve two targets with an angular separation of less than 25 °. The WSCMF method provided by the invention can distinguish two targets with an angle interval not less than 20 degrees. Although the root mean square error increases when the azimuth angle interval is less than 25 °, the estimation accuracy of the method of the present invention is the highest compared to the IMLSE method, the asmusch method, the IAA method, and the MUSIC method. In addition, it is easy to see that when the azimuth angle interval is greater than 50 °, the azimuth estimation accuracy of the WSCMF method, the IMLSE method, the asmuscic method, the IAA method, and the MUSIC method provided by the present invention is within an acceptable range, but the azimuth estimation accuracy of the WSCMF method provided by the present invention is the highest.
It can be found from fig. 5 that the estimation performance of the MUSIC method improves less as the signal-to-noise ratio increases. Although the performance of the IAA method, the IMLSE method and the asmusch method is improved, the azimuth estimation accuracy is not high enough when the signal-to-noise ratio is low. The WSCMF method provided by the invention adopts a regularization weighting sparse covariance matrix fitting method to estimate the noise power and the signal power, and well balances the relationship between the sparsity of a signal in a space domain and a fitting error through user parameters, so that the orientation estimation precision of the acoustic vector sensor array is improved when the signal-to-noise ratio is low.
As can be seen from fig. 6, the azimuth estimation performance of the MUSIC method and the IAA method, which do not consider the non-uniform noise, is seriously deteriorated as the worst noise power ratio increases. Although the root mean square error of the IMLSE method, the ASMUSIC method and the method provided by the invention is less affected with the increase of the worst noise power ratio, the root mean square error of the method provided by the invention is the smallest compared with the IMLSE method and the ASMUSIC method, which shows that the method provided by the invention has better azimuth estimation performance in a non-uniform noise scene.
While embodiments of the present invention have been described above, the above description is intended to be exemplary, not exhaustive, and not limited to the disclosed embodiments. Many modifications and variations will be apparent to those of ordinary skill in the art without departing from the scope and spirit of the illustrated embodiments.

Claims (4)

1. The method for estimating the direction of arrival of the acoustic vector sensor array under the nonuniform noise is characterized by comprising the following steps of:
step 1, receiving an output vector X output by an acoustic vector sensor array;
step 2, based on a sparse signal model, calculating a signal covariance matrix P and a noise covariance matrix Q according to the output vector X, and constructing a sparse covariance matrix R by using the obtained signal covariance matrix P and the noise covariance matrix Q;
step 3, constructing a cost function based on the sparse covariance matrix R, estimating sparse signal power by using the cost function, and obtaining a signal power vector
Figure FDA0003581551080000011
Step 4, for signal power vector
Figure FDA0003581551080000012
And searching a spectral peak, wherein the sound source position corresponding to the obtained spectral peak is the estimated target azimuth.
2. The method according to claim 1, wherein the step 2 comprises the following sub-steps:
Step 21, constructing a sparse signal model of the acoustic vector sensor array under the condition of non-uniform noise:
supposing that K far-field equipower narrowband signals with the wavelength of lambda are incident on a uniform linear array consisting of M acoustic vector sensors, the spacing between array elements is d, and the vibration velocity components of the narrowband signals along the directions of an x axis and a y axis are u respectivelyxAnd uyAnd the X-axis and the y-axis have orthogonality, the output vector X of the acoustic vector sensor array can be expressed as L in the case of fast beat number
X=A(θ)S+W
Wherein X is [ X (1), X (2), …, X (L) ]]Represents an output vector of the acoustic vector sensor array, and a (θ) ═ a (θ)1),a(θ2),…,a(θK)]An array manifold matrix representing an array of ideal acoustic vector sensors, θ ═ θ [ θ ]1,θ2,…,θK]The azimuth angle of the target is represented,
Figure FDA0003581551080000021
Figure FDA0003581551080000022
represents the kth orderThe target manifold vector is then used to determine,
Figure FDA0003581551080000023
representing the Cartesian product, h (θ)k)=[1 cosθk sinθk]TAn array manifold vector representing the incidence of the kth target on the acoustic vector transducer,
Figure FDA0003581551080000024
Figure FDA0003581551080000025
denotes a sound pressure response coefficient of the kth target incident on the mth acoustic vector sensor, W ═ W (1), W (2), …, W (l)]Represents a noise vector, S ═ S (1), S (2), …, S (l)]Representing a vector of input signal waveforms;
assuming that the angle space (-180 deg.) is uniformly divided into N discrete grids, and the number N of the discrete grids is far larger than the number K of the sound sources, the incoming wave directions of all the sound sources can be expressed as
Figure FDA0003581551080000026
The sparse signal model of the acoustic vector sensor array may be represented as
Figure FDA0003581551080000027
In the formula,
Figure FDA0003581551080000028
an array manifold matrix representing an array of acoustic vector sensors under a sparse signal model,
Figure FDA0003581551080000029
is a row sparse matrix and its non-zero rows represent the locations of the real signals;
step 22, calculating a signal covariance matrix P and a noise covariance matrix Q:
the signal covariance matrix P can be expressed as
Figure FDA00035815510800000210
Wherein E {. cndot } represents the desired operation, Diag {. cndot } represents the diagonal matrix operation,
Figure FDA00035815510800000211
represents the power of the nth grid corresponding signal;
the noise covariance matrix Q can be expressed as
Q=diag{q1,…,qM}
In the formula,
Figure FDA00035815510800000212
step 23, constructing a sparse covariance matrix R:
array manifold matrix defining a new acoustic vector transducer array
Figure FDA0003581551080000031
Wherein I3MRepresenting a 3M x 3M identity matrix, the sparse covariance matrix R can be represented as
R=DΞDH
Wherein,
Figure FDA0003581551080000032
in the formula, the first N terms of the diagonal elements of xi represent the signal power, and the N +1 to N +3M represent the noise power output by each channel in the acoustic vector sensor array. Thus, the noise covariance matrix Q can be further expressed as
Figure FDA0003581551080000033
In the formula, DnLine n data representing D, DN+l=ul,l=1,...,3M,ulIs represented by3MColumn l.
3. The method according to claim 2, wherein step 3 comprises the following sub-steps:
Step 31, constructing an interference-plus-noise covariance matrix B according to the sparse covariance matrix Rn
According to the sparse covariance matrix R constructed in step 23, the signal covariance matrix corresponding to the nth grid can be expressed as
Figure FDA0003581551080000041
Then the interference plus noise covariance matrix corresponding to the signal on the nth grid is
Bn=R-Rn=DΞDH-Rn
Step 32, based on the fitting criterion of the sparse covariance matrix, constructing the following cost function:
Figure FDA0003581551080000042
in the formula, | \ | non-countingFThe operation of taking the F norm is shown,
Figure FDA0003581551080000043
representing the sample covariance matrix, q representing the user parameter for the compensation signal, lambdasRegularization parameter representing the relationship of sparsity of the equilibrium signal to fitting error, ands>0;
step 33, solving the sparse signal power by using an iterative method, and obtaining a signal power vector
Figure FDA0003581551080000044
4. The method according to claim 3, wherein step 33 comprises the following sub-steps:
step 331, singular value decomposition is carried out on the sample covariance matrix to obtain
Figure FDA0003581551080000045
In the formula,
Figure FDA0003581551080000046
a signal sub-space is represented that is,
Figure FDA0003581551080000047
representing a noise subspace;
the cost function can be expressed as
Figure FDA0003581551080000048
Step 332, solving the sparse signal power by using an iterative method, setting the j +1 th iteration to obtain a unique solution, and obtaining the sparse signal power
Figure FDA0003581551080000051
In the formula, anIs composed of
Figure FDA0003581551080000052
In the form of a short hand of (a),
Figure FDA0003581551080000053
Figure FDA0003581551080000054
Figure FDA0003581551080000055
step 333, obtaining the signal power vector of the signals corresponding to the N grids according to the sparse signal power obtained in step 332 as
Figure FDA0003581551080000056
CN202210353033.9A 2022-04-06 2022-04-06 Method for estimating direction of arrival of acoustic vector sensor array under non-uniform noise Pending CN114755628A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210353033.9A CN114755628A (en) 2022-04-06 2022-04-06 Method for estimating direction of arrival of acoustic vector sensor array under non-uniform noise

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210353033.9A CN114755628A (en) 2022-04-06 2022-04-06 Method for estimating direction of arrival of acoustic vector sensor array under non-uniform noise

Publications (1)

Publication Number Publication Date
CN114755628A true CN114755628A (en) 2022-07-15

Family

ID=82328220

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210353033.9A Pending CN114755628A (en) 2022-04-06 2022-04-06 Method for estimating direction of arrival of acoustic vector sensor array under non-uniform noise

Country Status (1)

Country Link
CN (1) CN114755628A (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114966640A (en) * 2022-07-29 2022-08-30 宁波博海深衡科技有限公司 Direction estimation method and system based on array background noise statistical covariance estimation

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114966640A (en) * 2022-07-29 2022-08-30 宁波博海深衡科技有限公司 Direction estimation method and system based on array background noise statistical covariance estimation
CN114966640B (en) * 2022-07-29 2022-12-09 宁波博海深衡科技有限公司 Orientation estimation method and system based on array background noise statistics covariance estimation

Similar Documents

Publication Publication Date Title
CN108375763B (en) Frequency division positioning method applied to multi-sound-source environment
CN110113085B (en) Wave beam forming method and system based on covariance matrix reconstruction
CN107315162B (en) Far-field coherent signal DOA estimation method based on interpolation transformation and beam forming
CN110045321B (en) Robust DOA estimation method based on sparse and low-rank recovery
CN110045323B (en) Matrix filling-based co-prime matrix robust adaptive beamforming algorithm
CN109298383A (en) A kind of relatively prime battle array direction of arrival angle estimation method based on variational Bayesian
CN109557504B (en) Method for positioning near-field narrow-band signal source
CN108828586B (en) Bistatic MIMO radar angle measurement optimization method based on beam domain
CN108761394A (en) A kind of high-resolution low sidelobe based on space-time processing deconvolutes Power estimation method
CN108761380B (en) Target direction of arrival estimation method for improving precision
CN112180339A (en) Radar echo signal accurate direction finding method based on sparse processing
CN114814830A (en) Meter-wave radar low elevation height measurement method based on robust principal component analysis noise reduction
CN108594165B (en) Narrow-band signal direction-of-arrival estimation method based on expectation maximization algorithm
CN114755628A (en) Method for estimating direction of arrival of acoustic vector sensor array under non-uniform noise
Zuo et al. Improved capon estimator for high-resolution doa estimation and its statistical analysis
CN109541572B (en) Subspace orientation estimation method based on linear environment noise model
CN113093098B (en) Axial inconsistent vector hydrophone array direction finding method based on lp norm compensation
CN115980721A (en) Array self-correcting method for error-free covariance matrix separation
CN114167346B (en) DOA estimation method and system based on covariance matrix fitting array element expansion
CN116226611A (en) Chirp signal direction-of-arrival estimation method based on fractional domain deconvolution beam forming
CN114035157B (en) Sub-band delay estimation method and system based on expectation maximization algorithm
CN111323750B (en) Direct positioning method based on acoustic vector array network
CN109298381A (en) A kind of relatively prime battle array coherent signal azimuth estimation method based on variational Bayesian
CN114184999B (en) Method for processing generated model of cross-coupling small-aperture array
CN111431575B (en) Incoming wave direction sparse reconstruction method based on conventional beam forming

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