CN112731283A - High subsonic speed flying target acoustic direction finding method based on multistage wiener filter - Google Patents

High subsonic speed flying target acoustic direction finding method based on multistage wiener filter Download PDF

Info

Publication number
CN112731283A
CN112731283A CN202011553627.1A CN202011553627A CN112731283A CN 112731283 A CN112731283 A CN 112731283A CN 202011553627 A CN202011553627 A CN 202011553627A CN 112731283 A CN112731283 A CN 112731283A
Authority
CN
China
Prior art keywords
signal
axis
target
matrix
acoustic
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.)
Granted
Application number
CN202011553627.1A
Other languages
Chinese (zh)
Other versions
CN112731283B (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.)
Unite 91550 Of Pla
Original Assignee
Unite 91550 Of Pla
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 Unite 91550 Of Pla filed Critical Unite 91550 Of Pla
Priority to CN202011553627.1A priority Critical patent/CN112731283B/en
Publication of CN112731283A publication Critical patent/CN112731283A/en
Application granted granted Critical
Publication of CN112731283B publication Critical patent/CN112731283B/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/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/86—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 with means for eliminating undesired waves, e.g. disturbing noises
    • 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
    • G01S3/8027—By vectorial composition of signals received by plural, differently-oriented transducers
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Abstract

The invention discloses a high subsonic speed flying target acoustic direction finding method based on a multistage wiener filter, which adopts an acoustic sensor array uniformly distributed on a spatial cross array to carry out direction finding, firstly, a broadband time domain array signal model is established to obtain a received signal matrix; performing wavelet transformation on each discrete received signal to obtain a wavelet matrix of the received signal matrix under a certain scale; performing multi-stage wiener filtering on the wavelet-transformed signal by using a multi-stage wiener filter to obtain a signal subspace and a noise subspace, constructing a two-dimensional spatial spectrum function, and combining two variables of a pitch angle and an azimuth angle
Figure DDA0002858677470000011
And according to the value range and the set step length, carrying out two-dimensional spectral peak search on the two-dimensional spatial spectral function of the target sound signal to obtain the optimal estimation value of the two-dimensional arrival angle of the target sound signal. Hair brushThe method effectively reduces the calculation complexity of direction finding, improves the real-time performance of high subsonic flight target estimation, and expands the application range of the DOA estimation method based on the multistage wiener filter.

Description

High subsonic speed flying target acoustic direction finding method based on multistage wiener filter
Technical Field
The invention relates to the field of target detection and tracking, in particular to a high subsonic speed flying target acoustic direction finding method based on a multistage wiener filter.
Background
The over-the-horizon detection of the sea surface low-altitude flying target is difficult in the conventional radar detection technology due to the restriction of the height of the sea surface platform. Meanwhile, due to the adverse factors of complex electromagnetic environment near the sea surface platform, serious low elevation sea clutter, multipath effect influence and the like, the radar detection of the sea surface short-range low-altitude flying target also has certain problems. The acoustic detection method detects the sea surface targets by receiving the flight noise of the sea surface targets, has the advantages of being free from the influence of electromagnetic and light environments, omni-directional monitoring, free from the restriction of sight conditions and the like, and particularly has unique background and environmental noise and less clutter interference when being applied to a sea surface platform, so that the acoustic detection method has unique advantages in the aspect of detecting the sea surface low-altitude flight targets.
With the sound velocity as a boundary, low-altitude flight targets can be divided into two categories, namely supersonic velocity targets and hypersonic velocity targets. The supersonic velocity target can generate disturbance to the surrounding air to form Mach waves when flying, the wave front of the Mach waves is called shock waves, the shock waves are propagated in a conical surface mode, and the time domain characteristics are obvious, so that relatively mature theories and technologies are provided for acoustic tracking of the supersonic velocity target. For a high subsonic speed flying target, the high subsonic speed flying target presents broadband, transient and non-stationary characteristics, and the acoustic tracking method for the high subsonic speed flying target is relatively less researched.
For the acoustic detection method, because it does not have an active detection capability for a target distance, and a single stationary acoustic array has been proved to be unable to estimate a target motion trajectory only by means of a Direction of arrival, and the target motion parameter estimation needs acoustic sensors in the array to reach a distance of hundreds of meters, the Direction of arrival (DOA) estimation for target noise is mainly realized by a single sensor at present, and the target trajectory estimation is realized by means of angle intersection and the like on the basis. Meanwhile, in a sea surface low-altitude flying target detection scene, the acoustic detection method is limited by the sound wave frequency and wavelength, and the DOA estimation accuracy is also limited, so that the target azimuth information obtained by the acoustic detection method is mainly used for providing guide information for other high-accuracy detection systems, and higher requirements are provided for the real-time performance of the acoustic DOA estimation method.
From the current research situation of the acoustic direction finding method of the current broadband signal, the currently proposed method comprises a broadband signal subspace method, a broadband beam domain method, a broadband cyclostationary DOA estimation method and a broadband array signal DOA estimation based on a time domain model. The broadband signal subspace method needs to invert the signal correlation matrix, thereby increasing the calculation complexity; the broadband wave beam domain method mainly aims at scenes in which multiple correlation signals are incident simultaneously, the broadband cyclostationary DOA estimation method generally aims at continuously occurring communication signals, and sea surface high subsonic speed flying target noise signals are short in duration and do not have the cyclostationary characteristic. The DOA estimation of the broadband array signal based on the time domain model can be simultaneously suitable for narrowband and broadband signals, a Markov chain Monte Carlo method is generally adopted for solving, and the problem of high calculation complexity also exists.
Disclosure of Invention
The invention discloses a method for acoustically direction finding of a high subsonic speed flying target, which aims at solving the problems of low real-time performance and limited precision of the conventional acoustic direction finding of the high subsonic speed flying target and is characterized in that a sound sensor array uniformly distributed on a spatial cross array is adopted for direction finding, the spatial cross array comprises three installation shafts of an X axis, a Y axis and a Z axis which are mutually vertical, the intersection point of the three installation shafts is called an O point, a sound sensor positioned at the O point is called a reference sensor, and 2M non-reference sensors are uniformly arranged on each installation shaft except the reference sensor;
when the direction finding is carried out on the high subsonic speed flying target, the method specifically comprises the following steps:
s1, aiming at the acoustic sensor array, establishing a broadband time domain array signal model to obtain a received signal matrix;
s2, performing wavelet transformation on each discrete received signal to obtain a wavelet matrix of the received signal matrix under a certain scale;
s3, according to the wavelet transformation result of step S2, selecting the wavelet transformation result of the discrete signal vector S (n) received by the reference sensor as the reference signal d0(k) Selecting a wavelet matrix WT (a) obtained by wavelet transforming a received signal matrix0B) is the initial received signal, a0For wavelet transformation of scale, multilevel wiener filtering is performed on the two signals by using a multilevel wiener filter to obtain a signal subspace USSum noise subspace UNK is the serial number of the sampling point, and n is the moment of the sampling point;
s4, constructing a two-dimensional spatial spectrum function by using the signal subspace or the noise subspace
Figure BDA0002858677450000021
For the combination of two variables of pitch angle and azimuth angle
Figure BDA0002858677450000022
According to the value range and the set step length, two-dimensional spectral peak search is carried out on the two-dimensional spatial spectral function of the target sound signal, and the optimal estimated value of the two-dimensional angle of arrival of the target sound signal is obtained as follows:
Figure BDA0002858677450000031
two-dimensional angle of arrival optimal estimation value
Figure BDA0002858677450000032
The direction finding result of the high subsonic speed flying target is obtained, and therefore the direction finding of the high subsonic speed flying target is completed.
When the method is used for carrying out direction finding on the high subsonic speed flying target, the specific steps comprise:
s1, aiming at the acoustic sensor array, establishing a broadband time domain array signal model to obtain a received signal matrix;
for a broadband incident signal s (t), its frequencyIn the range of [ fl,fh]With a bandwidth of B, fh＝fl+ B, in order to avoid the direction finding blur, the distance d between the acoustic sensors should satisfy:
d≤c/(2fh),
where c is the speed of sound in air.
The broadband target sound source signal received by the reference sensor is s (t), and for the signal s (t-tau) received by the non-reference sensor, a band-pass signal reconstruction and sampling method is adopted, and the method is expressed as follows:
Figure BDA0002858677450000033
where τ is the time delay of the non-reference sensor received signal relative to the reference sensor received signal, TsFor the sampling period ψ (t) is the band pass sampled signal, which takes the band pass form of a sinc function or a band pass ellipsoid wave function, L denotes the length of the window function and L denotes the delay of the window function.
And (3) performing time domain discrete sampling on s (t-tau), and replacing l with n-l to obtain the output of the sensor receiving signal at the n moment as follows:
Figure BDA0002858677450000034
wherein the bandpass window function hlAnd (τ) ═ ψ (l- τ) ω (l- τ), ω (l- τ) is a window function.
For 2M acoustic sensors on a certain installation axis, the received signals can be expressed as:
Figure BDA0002858677450000035
where m is the number of delays. In case n and m take different values, the received signals of different array elements are represented as a received signal matrix Sτ(n)：
Sτ(n)＝H(n,τ)s(n)，
Wherein H (n, tau) is a time delay difference matrixS (n) is a discrete signal vector obtained by sampling the signal received by the reference sensor and intercepting the signal at L-1 sampling moments respectively expanded to the front and the back by taking the n moment as the center, and receiving a signal matrix SτEach row in (n) represents a discrete received signal of one acoustic sensor. For matrix H (n, τ), the expression is:
Figure BDA0002858677450000041
and S2, performing wavelet transformation on each discrete received signal to obtain a wavelet matrix of the received signal matrix under a certain scale.
Specifically, for the received signal matrix S in step S1τ(n) performing wavelet transform to obtain a received signal matrix Sτ(n) in a0Wavelet matrix WT (a) at scale0B), the expression is:
WT(a0,b)＝A(a0,θ)Y(a0,b)，
wherein A (a)0τ) is the result of a wavelet transform of matrix H (n, τ), Y (a)0B) wavelet transform results of discrete signal vectors s (n);
s3, according to the wavelet transformation result of step S2, selecting the wavelet transformation result of the discrete signal vector S (n) received by the reference sensor as the reference signal d0(k) Selecting a wavelet matrix WT (a) obtained by wavelet transforming a received signal matrix0And b) initially receiving the signals, and performing multistage wiener filtering on the two signals by using a multistage wiener filter to obtain a signal subspace USK is the sampling point serial number; the specific process of the multistage wiener filtering comprises the calculation of a reference signal d0(k) Cross correlation vector with initial received signal
Figure BDA0002858677450000042
Order to
Figure BDA0002858677450000043
Let i be 1,2,3,.., 6M, and loop forward recursion:
Figure BDA0002858677450000044
Figure BDA0002858677450000045
Xi(k)＝Xi-1(k)-hidi(k)，
wherein h isiFor the sparse vector obtained for the ith recursion, di(k) For the reference signal, X, obtained in the ith recursioni(k) Receiving signals obtained for the ith recursion;
computing signal subspace USSum noise subspace UN：
US＝span{h1,h2,…,hP+2}，
UN＝span{hp+3,hp+4,…,h6M}，
Wherein p is the number of high subsonic flight targets;
s4, constructing a two-dimensional spatial spectrum function by using the signal subspace or the noise subspace
Figure BDA0002858677450000051
Specifically, the two-dimensional spatial spectrum function constructed using the noise subspace is:
Figure BDA0002858677450000052
the two-dimensional spatial spectrum function constructed using the signal subspace is:
Figure BDA0002858677450000053
wherein the content of the first and second substances,
Figure BDA0002858677450000054
for target acoustic signalsAnd the value of the direction vector is determined by the relative position of each non-reference acoustic sensor and the reference sensor.
For the combination of two variables of pitch angle and azimuth angle
Figure BDA0002858677450000055
According to the value range and the set step length, two-dimensional spectral peak search is carried out on the two-dimensional spatial spectral function of the target sound signal, and the optimal estimated value of the two-dimensional angle of arrival of the target sound signal is obtained as follows:
Figure BDA0002858677450000056
two-dimensional angle of arrival optimal estimation value
Figure BDA0002858677450000057
The direction finding result of the high subsonic speed flying target is obtained, and therefore the direction finding of the high subsonic speed flying target is completed.
The space cross array comprises three installation axles of mutually perpendicular's X axle, Y axle and Z axle, and three installation axles of X axle, Y axle and Z axle are space rectangular coordinate and distribute, X axle and Y axle are located the horizontal plane, the Z axle is perpendicular with the horizontal plane, Z axle forward directional sky. When the azimuth angle and the pitch angle of the incident wave of the high subsonic flight target are respectively
Figure BDA0002858677450000058
And when theta is satisfied, the pitch angle theta is an included angle between a connecting line from the reference sensor to the target and the positive direction of the Z axis, the transmission delay tau of the received signals between the adjacent sensors on the same installation shaft is alpha d/c, alpha is a delay factor, and the values of the corresponding alpha are respectively the values of the X axis, the Y axis and the Z axis
Figure BDA0002858677450000061
And cos θ, direction vector for target acoustic signal
Figure BDA0002858677450000062
The expression is as follows:
Figure BDA0002858677450000063
wherein, ax、ayAnd azTarget acoustic signal direction vectors of an X axis, a Y axis and a Z axis respectively have the following expressions:
Figure BDA0002858677450000064
the invention has the beneficial effects that:
(1) according to the method, by constructing a time domain broadband signal model of a uniform linear space cross array and adopting a time domain broadband DOA (direction of arrival) rapid estimation method based on MSWF (minimum shift function), the calculation complexity of the system is effectively reduced, and the real-time performance of the estimation of the high subsonic speed flight target is improved.
(2) The method effectively expands the application range of the DOA estimation method based on the multistage wiener filter by constructing the broadband receiving model in the time domain, and provides a technical basis for effectively realizing the acoustic detection of the offshore high subsonic flight target.
(3) The method realizes effective suppression of environmental noise by performing wavelet transformation on the received signals, and meanwhile, gathers limited time domain signals, thereby improving the anti-noise performance and direction finding precision of the method.
Drawings
FIG. 1 is a schematic view of a uniform linear spatial cross array.
FIG. 2 shows the comparison result of the simulation of the actual calculated quantity of the method of the present invention and the MUSIC method.
Detailed Description
For a better understanding of the present disclosure, an example is given here.
The embodiment discloses a high subsonic flight target acoustic direction finding method based on a Multi-Stage Wiener Filtering (MSWF). the direction finding is carried out by adopting an acoustic sensor array uniformly distributed on a spatial cross array, the spatial cross array comprises three installation shafts of an X axis, a Y axis and a Z axis which are mutually perpendicular, the intersection point of the three installation shafts is called an O point, an acoustic sensor positioned at the O point is called a reference sensor, and 2M non-reference sensors are uniformly installed on each installation shaft except the reference sensor.
When the direction finding is carried out on the high subsonic speed flying target, the method specifically comprises the following steps:
s1, aiming at the acoustic sensor array, establishing a broadband time domain array signal model to obtain a received signal matrix;
for a broadband incident signal s (t) with a frequency range of [ f ]l,fh]With a bandwidth of B, fh＝fl+ B, in order to avoid the direction finding blur, the distance d between the acoustic sensors should satisfy:
d≤c/(2fh),
where c is the speed of sound in air.
The broadband target sound source signal received by the reference sensor is s (t), and for the signal s (t-tau) received by the non-reference sensor, a band-pass signal reconstruction and sampling method is adopted, and the method is expressed as follows:
Figure BDA0002858677450000071
where τ is the time delay of the non-reference sensor received signal relative to the reference sensor received signal, TsFor the sampling period ψ (t) is the band pass sampled signal, which takes the band pass form of a sinc function or a band pass ellipsoid wave function, L denotes the length of the window function and L denotes the delay of the window function.
And (3) performing time domain discrete sampling on s (t-tau), and replacing l with n-l to obtain the output of the sensor receiving signal at the n moment as follows:
Figure BDA0002858677450000072
wherein the bandpass window function hlAnd (τ) ═ ψ (l- τ) ω (l- τ), ω (l- τ) is a window function.
For 2M acoustic sensors on a certain installation axis, the received signals can be expressed as:
Figure BDA0002858677450000073
where m is the number of delays. In case n and m take different values, the received signals of different array elements are represented as a received signal matrix Sτ(n)：
Sτ(n)＝H(n,τ)s(n)，
Wherein, H (n, tau) is a time delay difference matrix, S (n) is a discrete signal vector obtained by sampling the signal received by the reference sensor, intercepting the signal by expanding L-1 sampling moments forward and backward respectively by taking n moment as the center, and receiving a signal matrix SτEach row in (n) represents a discrete received signal of one acoustic sensor. For matrix H (n, τ), the expression is:
Figure BDA0002858677450000081
and S2, performing wavelet transformation on each discrete received signal to obtain a wavelet matrix of the received signal matrix under a certain scale.
Specifically, for the received signal matrix S in step S1τ(n) performing wavelet transform to obtain a received signal matrix Sτ(n) in a0Wavelet matrix WT (a) at scale0B), the expression is:
WT(a0,b)＝A(a0,θ)Y(a0,b)，
wherein A (a)0τ) is the result of a wavelet transform of matrix H (n, τ), Y (a)0B) wavelet transform results of discrete signal vectors s (n);
s3, according to the wavelet transformation result of step S2, selecting the wavelet transformation result of the discrete signal vector S (n) received by the reference sensor as the reference signal d0(k) Selecting a wavelet matrix WT (a) obtained by wavelet transforming a received signal matrix0And b) initially receiving the signals, and performing multistage wiener filtering on the two signals by using a multistage wiener filter to obtain a signal subspace USK is the sampling point serial number; the specific process of the multistage wiener filtering comprises the calculation of a reference signal d0(k) Cross correlation vector with initial received signal
Figure BDA0002858677450000082
Order to
Figure BDA0002858677450000083
Let i be 1,2,3,.., 6M, and loop forward recursion:
Figure BDA0002858677450000084
Figure BDA0002858677450000085
Xi(k)＝Xi-1(k)-hidi(k)，
wherein h isiFor the sparse vector obtained for the ith recursion, di(k) For the reference signal, X, obtained in the ith recursioni(k) Receiving signals obtained for the ith recursion;
computing signal subspace USSum noise subspace UN：
US＝span{h1,h2,…,hP+2}，
UN＝span{hp+3,hp+4,…,h6M}，
Wherein p is the number of high subsonic flight targets;
s4, constructing a two-dimensional spatial spectrum function by using the signal subspace or the noise subspace
Figure BDA0002858677450000091
Specifically, the two-dimensional spatial spectrum function constructed using the noise subspace is:
Figure BDA0002858677450000092
the two-dimensional spatial spectrum function constructed using the signal subspace is:
Figure BDA0002858677450000093
wherein the content of the first and second substances,
Figure BDA0002858677450000094
the direction vector of the target acoustic signal is determined by the relative position of each non-reference acoustic sensor and the reference sensor.
For the combination of two variables of pitch angle and azimuth angle
Figure BDA0002858677450000095
According to the value range and the set step length, two-dimensional spectral peak search is carried out on the two-dimensional spatial spectral function of the target sound signal, and the optimal estimated value of the two-dimensional angle of arrival of the target sound signal is obtained as follows:
Figure BDA0002858677450000096
two-dimensional angle of arrival optimal estimation value
Figure BDA0002858677450000097
The direction finding result of the high subsonic speed flying target is obtained, and therefore the direction finding of the high subsonic speed flying target is completed.
The space cross array is composed of three installation shafts of an X shaft, a Y shaft and a Z shaft which are perpendicular to each other, the three installation shafts of the X shaft, the Y shaft and the Z shaft are distributed in a space rectangular coordinate mode, the X shaft and the Y shaft are located on a horizontal plane, the Z shaft is perpendicular to the horizontal plane, the forward direction of the Z shaft points to the sky, the drawing 1 is a schematic diagram of the uniform linear space cross array, 101 is a high subsonic speed flight target, and 102 is the space cross array. When the azimuth angle and the pitch angle of the incident wave of the high subsonic flight target are respectively
Figure BDA0002858677450000098
And theta, the pitch angle theta isThe included angle between the connecting line from the reference sensor to the target and the Z axis in the forward direction is that for the transmission delay tau between the adjacent sensors on the same installation shaft, alpha is a delay factor, and for the X axis, the Y axis and the Z axis, the corresponding alpha values are respectively
Figure BDA0002858677450000101
And cos θ, direction vector for target acoustic signal
Figure BDA0002858677450000102
The expression is as follows:
Figure BDA0002858677450000103
wherein, ax、ayAnd azTarget acoustic signal direction vectors of an X axis, a Y axis and a Z axis respectively have the following expressions:
Figure BDA0002858677450000104
under the conditions that the number of array elements is 4, 8, 12, 16 and 20, the target is 1, the signal-to-noise ratio is 0dB, the snapshot number is 100 and 200, the spacing of the array elements is 0.5, and the search step length is 0.5 degrees, the simulation is carried out on the machine period number required by the operation of the method and the classical MUSIC method. FIG. 2 shows the comparison result of the actual computation simulation of the method of the present invention and the classical MUSIC method. When the number of array elements is less, considering relevant processes such as variable definition, storage and the like, the actual machine cycle number required by the two methods is basically the same; the actual number of machine cycles required for the method of the invention is less as the number of array elements increases. At the same time, as the number of fast beats decreases, the actual number of machine cycles required for both methods also gradually decreases.
The above description is only an example of the present application and is not intended to limit the present application. Various modifications and changes may occur to those skilled in the art. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present application should be included in the scope of the claims of the present application.

Claims (7)

1. A high subsonic speed flying target acoustic direction finding method based on a multistage wiener filter is characterized in that a sound sensor array uniformly distributed on a spatial cross array is adopted for direction finding, the spatial cross array comprises three installation shafts of an X axis, a Y axis and a Z axis which are mutually perpendicular, the intersection point of the three installation shafts is called an O point, a sound sensor positioned at the O point is called a reference sensor, and 2M non-reference sensors are uniformly installed on each installation shaft except the reference sensor;
when the direction finding is carried out on the high subsonic speed flying target, the method specifically comprises the following steps:
s1, aiming at the acoustic sensor array, establishing a broadband time domain array signal model to obtain a received signal matrix;
s2, performing wavelet transformation on each discrete received signal to obtain a wavelet matrix of the received signal matrix under a certain scale;
s3, according to the wavelet transformation result of step S2, selecting the wavelet transformation result of the discrete signal vector S (n) received by the reference sensor as the reference signal d0(k) Selecting a wavelet matrix WT (a) obtained by wavelet transforming a received signal matrix0B) is the initial received signal, a0For wavelet transformation of scale, multilevel wiener filtering is performed on the two signals by using a multilevel wiener filter to obtain a signal subspace USSum noise subspace UNK is the serial number of the sampling point, and n is the moment of the sampling point;
s4, constructing a two-dimensional spatial spectrum function by using the signal subspace or the noise subspace
Figure FDA0002858677440000011
For the combination of two variables of pitch angle and azimuth angle
Figure FDA0002858677440000012
According to the value range and the set step length, two-dimensional spectral peak search is carried out on the two-dimensional space spectral function of the target sound signal to obtain the optimal estimation of the two-dimensional angle of arrival of the target sound signalThe values are:
Figure FDA0002858677440000013
two-dimensional angle of arrival optimal estimation value
Figure FDA0002858677440000014
The direction finding result of the high subsonic speed flying target is obtained, and therefore the direction finding of the high subsonic speed flying target is completed.
2. The method for acoustically measuring the direction of a target flying at high subsonic speed based on multistage wiener filters according to claim 1, wherein the frequency range of the broadband incident signal s (t) is [ f [ ]l,fh]With a bandwidth of B, fh＝fl+ B, in order to avoid the direction finding fuzzy, the acoustic sensor interval d on the same installation shaft should satisfy:
d≤c/(2fh),
where c is the speed of sound in air.
3. The acoustic direction finding method for the flying target with high subsonic speed based on the multistage wiener filter as claimed in claim 1 or 2, wherein said step S1 specifically includes that the reference sensor receives a broadband target sound source signal S (t), and for a non-reference sensor received signal S (t- τ), a band-pass signal reconstruction and sampling method is adopted, which is expressed as:
Figure FDA0002858677440000021
where τ is the time delay of the non-reference sensor received signal relative to the reference sensor received signal, TsFor the sampling period, ψ (t) is a band-pass sampling signal which adopts a band-pass form of a sinc function or a band-pass elliptic spherical wave function, L represents the length of a window function, and L represents the delay of the window function;
and (3) performing time domain discrete sampling on s (t-tau), and replacing l with n-l to obtain the output of a sensor receiving signal at a sampling point n moment as follows:
Figure FDA0002858677440000022
wherein the bandpass window function hl(τ) ═ ψ (l- τ) ω (l- τ), ω (l- τ) is a window function;
for 2M acoustic sensors on a certain installation axis, the received signals can be expressed as:
Figure FDA0002858677440000023
wherein m is the number of delays; in case n and m take different values, the received signals of different array elements are represented as a received signal matrix Sτ(n)：
Sτ(n)＝H(n,τ)s(n)，
Wherein, H (n, tau) is a time delay difference matrix, S (n) is a time delay difference matrix, and after sampling the signal received by the reference sensor, the signal is intercepted by expanding L-1 sampling moments to the front and the back by taking n moment as a center to obtain a discrete signal vector, and a received signal matrix Sτ(n) each row represents a discrete received signal of one acoustic sensor; for matrix H (n, τ), the expression is:
Figure FDA0002858677440000031
4. the acoustic direction finding method for high subsonic flight targets based on multistage wiener filters as claimed in claim 3, wherein said step S2 is performed with respect to the received signal matrix S1τ(n) performing wavelet transform to obtain a received signal matrix Sτ(n) in a0Wavelet matrix WT (a) at scale0B), the expression is:
WT(a0,b)＝A(a0,θ)Y(a0,b)，
wherein A (a)0τ) is the result of a wavelet transform of matrix H (n, τ), Y (a)0And b) is the result of wavelet transform of a discrete signal vector s (n).
5. The acoustic direction finding method for the high subsonic flying target based on the multistage wiener filter as claimed in claim 3, wherein the multistage wiener filter of step S3 includes calculating the reference signal d0(k) Cross correlation vector with initial received signal
Figure FDA0002858677440000032
Order to
Figure FDA0002858677440000033
Let i be 1,2,3,.., 6M, and loop forward recursion:
Figure FDA0002858677440000034
Figure FDA0002858677440000035
Xi(k)＝Xi-1(k)-hidi(k)，
wherein h isiFor the sparse vector obtained for the ith recursion, di(k) For the reference signal, X, obtained in the ith recursioni(k) Receiving signals obtained for the ith recursion;
computing signal subspace USSum noise subspace UN：
US＝span{h1,h2,…,hP+2}，
UN＝span{hp+3,hp+4,…,h6M}，
Where p is the number of high subsonic flight targets.
6. The method for acoustically direction finding of a high subsonic velocity flying target based on multistage wiener filters according to claim 3, wherein said two-dimensional spatial spectral function is constructed
Figure FDA0002858677440000041
The two-dimensional spatial spectrum function constructed using the noise subspace is:
Figure FDA0002858677440000042
the two-dimensional spatial spectrum function constructed using the signal subspace is:
Figure FDA0002858677440000043
wherein the content of the first and second substances,
Figure FDA0002858677440000044
is the direction vector of the target acoustic signal,
Figure FDA0002858677440000045
the pitch angle and the azimuth angle of the high subsonic flight target are determined by the relative positions of the non-reference acoustic sensors and the reference sensor.
7. The acoustic direction finding method for the target flying at high subsonic speed based on the multistage wiener filter as claimed in claim 1 or 2, wherein the spatial cross array is composed of three installation shafts of X axis, Y axis and Z axis which are perpendicular to each other, the three installation shafts of X axis, Y axis and Z axis are distributed in a space rectangular coordinate mode, the X axis and Y axis are located on a horizontal plane, the Z axis is perpendicular to the horizontal plane, and the Z axis is directed towards the sky; when the azimuth angle and the pitch angle of the incident wave of the high subsonic flight target are respectively
Figure FDA00028586774400000410
And when theta is satisfied, the pitch angle theta is an included angle between a connecting line from the reference sensor to the target and the positive direction of the Z axis, the transmission delay tau of the received signals between the adjacent sensors on the same installation shaft is alpha d/c, alpha is a delay factor, and the values of the corresponding alpha are respectively the values of the X axis, the Y axis and the Z axis
Figure FDA0002858677440000046
And cos θ, direction vector for target acoustic signal
Figure FDA0002858677440000047
The expression is as follows:
Figure FDA0002858677440000048
wherein, ax、ayAnd azTarget acoustic signal direction vectors of an X axis, a Y axis and a Z axis respectively have the following expressions:
Figure FDA0002858677440000049
CN202011553627.1A 2020-12-24 2020-12-24 High subsonic flight target acoustic direction finding method based on multistage wiener filter Active CN112731283B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011553627.1A CN112731283B (en) 2020-12-24 2020-12-24 High subsonic flight target acoustic direction finding method based on multistage wiener filter

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011553627.1A CN112731283B (en) 2020-12-24 2020-12-24 High subsonic flight target acoustic direction finding method based on multistage wiener filter

Publications (2)

Publication Number Publication Date
CN112731283A true CN112731283A (en) 2021-04-30
CN112731283B CN112731283B (en) 2023-07-11

Family

ID=75615502

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011553627.1A Active CN112731283B (en) 2020-12-24 2020-12-24 High subsonic flight target acoustic direction finding method based on multistage wiener filter

Country Status (1)

Country Link
CN (1) CN112731283B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114757242A (en) * 2022-06-16 2022-07-15 中国空气动力研究与发展中心低速空气动力研究所 Helicopter noise enhancement method and detection method based on cyclic wiener filtering

Citations (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20020158801A1 (en) * 2001-04-27 2002-10-31 Crilly William J. Wireless packet switched communication systems and networks using adaptively steered antenna arrays
US20040178951A1 (en) * 2002-03-13 2004-09-16 Tony Ponsford System and method for spectral generation in radar
CN102353947A (en) * 2011-07-08 2012-02-15 哈尔滨工程大学 Method for estimating target echo signal subspaces of passive radars based on CSA-MWF (correlation subtraction algorithm-multistage wiener filter)
JP2013152112A (en) * 2012-01-24 2013-08-08 Mitsubishi Electric Corp Time difference orientation detection device
CN103969626A (en) * 2014-05-20 2014-08-06 西安电子科技大学 Wideband digital wave beam forming method based on all-pass type variable fractional delay filter
JP2014142279A (en) * 2013-01-24 2014-08-07 Furuno Electric Co Ltd Signal processor, radar device, underwater detection device, and signal processing method
CN105158751A (en) * 2015-08-29 2015-12-16 许昌学院 Acoustic vector array fast DOA (Direction of Arrival) estimation method
US20150362580A1 (en) * 2014-06-11 2015-12-17 Agency For Defense Development Method for estimating angle of arrival of multi-target moving at high speed
CN107255793A (en) * 2017-06-16 2017-10-17 中国电子科技集团公司第二十九研究所 A kind of array df method and device for wideband OFDM signal of communication
CN108983144A (en) * 2018-05-29 2018-12-11 西北工业大学 It improves Wiener filter and carries out the estimation method of target bearing based on the filter
CN111175693A (en) * 2020-01-19 2020-05-19 河北科技大学 Direction-of-arrival estimation method and direction-of-arrival estimation device

Patent Citations (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20020158801A1 (en) * 2001-04-27 2002-10-31 Crilly William J. Wireless packet switched communication systems and networks using adaptively steered antenna arrays
US20040178951A1 (en) * 2002-03-13 2004-09-16 Tony Ponsford System and method for spectral generation in radar
CN102353947A (en) * 2011-07-08 2012-02-15 哈尔滨工程大学 Method for estimating target echo signal subspaces of passive radars based on CSA-MWF (correlation subtraction algorithm-multistage wiener filter)
JP2013152112A (en) * 2012-01-24 2013-08-08 Mitsubishi Electric Corp Time difference orientation detection device
JP2014142279A (en) * 2013-01-24 2014-08-07 Furuno Electric Co Ltd Signal processor, radar device, underwater detection device, and signal processing method
CN103969626A (en) * 2014-05-20 2014-08-06 西安电子科技大学 Wideband digital wave beam forming method based on all-pass type variable fractional delay filter
US20150362580A1 (en) * 2014-06-11 2015-12-17 Agency For Defense Development Method for estimating angle of arrival of multi-target moving at high speed
CN105158751A (en) * 2015-08-29 2015-12-16 许昌学院 Acoustic vector array fast DOA (Direction of Arrival) estimation method
CN107255793A (en) * 2017-06-16 2017-10-17 中国电子科技集团公司第二十九研究所 A kind of array df method and device for wideband OFDM signal of communication
CN108983144A (en) * 2018-05-29 2018-12-11 西北工业大学 It improves Wiener filter and carries out the estimation method of target bearing based on the filter
CN111175693A (en) * 2020-01-19 2020-05-19 河北科技大学 Direction-of-arrival estimation method and direction-of-arrival estimation device

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
宋爱民等: "非圆信号多级维纳滤波DOA估计求根算法", 电子科技大学学报, no. 01 *
胡政权等: "复解析小波变换的圆阵波达方向估计方法", 探测与控制学报, no. 04 *
谢纪岭;司锡才;唐建红: "基于多级维纳滤波器的二维测向算法及DSP实现", 宇航学报, no. 01 *
马社祥;赵朝霞: "空间交叉阵被动声定位算法及其性能分析", 系统工程与电子技术, no. 03 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114757242A (en) * 2022-06-16 2022-07-15 中国空气动力研究与发展中心低速空气动力研究所 Helicopter noise enhancement method and detection method based on cyclic wiener filtering
CN114757242B (en) * 2022-06-16 2022-09-23 中国空气动力研究与发展中心低速空气动力研究所 Helicopter noise enhancement method and detection method based on cyclic wiener filtering

Also Published As

Publication number Publication date
CN112731283B (en) 2023-07-11

Similar Documents

Publication Publication Date Title
CN106842128B (en) The acoustics tracking and device of moving target
CN103760529A (en) Efficient cascading space-time adaptive processing method based on passive detection
CN110824439B (en) Radar target rapid long-time coherent accumulation method
CN112346030A (en) Super-resolution direction-of-arrival estimation method for unmanned aerial vehicle cluster
CN112834980A (en) Vortex electromagnetic wave super-resolution direction finding method based on propagation operator
CN112731283B (en) High subsonic flight target acoustic direction finding method based on multistage wiener filter
CN113671486A (en) Radar detection system and method for rotor wing target
CN115856813A (en) Radar target sidelobe suppression method based on APC and IARFT cascade processing
CN112731292B (en) Low-altitude flight target signal time delay estimation method with local IMF energy weighting
CN115436898A (en) High-speed target coherent accumulation detection method for airborne distributed MIMO radar
CN111610505A (en) Airborne bistatic radar STAP algorithm based on exponential form time-varying weighting
CN115291207A (en) Multi-target detection method for small rotor unmanned aerial vehicle based on MIMO radar
CN110231590B (en) Array target angle measurement method based on DFT (discrete Fourier transform)
Wan et al. A range-Doppler-angle estimation method for passive bistatic radar
CN113075633A (en) Target positioning method based on distance-pitching joint dictionary
CN113447919A (en) Extended Kalman prediction angle tracking method
CN113126029A (en) Multi-sensor pulse sound source positioning method suitable for deep sea reliable acoustic path environment
CN111796251A (en) Method for processing time division MIMO radar signal by using GPU
Zima et al. Investigation and modeling of the spatiotemporal signal recorded by the digital antenna array
CN113671477B (en) Radar target distance estimation method based on graph signal processing
Raia et al. GSP-based DoA estimation for a multimission radar
CN111007472B (en) Clutter echo modeling method for hypersonic platform in complex motion state
CN116990791B (en) Multi-receiving-array synthetic aperture sonar echo simulation method
CN114779198B (en) Conformal array airborne radar space-time clutter spectrum adaptive compensation and clutter suppression method
CN115267721B (en) Ground moving target radial velocity estimation method based on double-frequency SAR

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