CN109799495A - A kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING - Google Patents
A kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING Download PDFInfo
- Publication number
- CN109799495A CN109799495A CN201910001466.6A CN201910001466A CN109799495A CN 109799495 A CN109799495 A CN 109799495A CN 201910001466 A CN201910001466 A CN 201910001466A CN 109799495 A CN109799495 A CN 109799495A
- Authority
- CN
- China
- Prior art keywords
- array element
- spectrum
- frequency
- power spectrum
- signal
- 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
Links
Abstract
The invention discloses a kind of broadband delay time estimation methods for high-fidelity ARRAY PROCESSING, include the following steps: that (1) seeks beam energy figure, estimation target bearing (2) reconstructs source signal in advance, detects the frequency domain signal-to-noise ratio of the power spectrum respective frequencies line spectrum of power spectrum line spectrum position (3) the calculating reference array element acquisition data of pre- reconstruct source signal；(4) phase difference between the array element of calculating respective frequencies line spectrum；(5) fusion is weighted between phase difference the array element estimated based on different frequency line spectrum, estimation echo signal broadband time delay vector (6) carries out secondary reconstruct to source signal from basic matrix acquisition data.Time delay between line spectrum of this method based on echo signal different frequency array element having the same, fusion is weighted to phase difference between the array element estimated based on different frequency line spectrum according to the frequency of line spectrum and frequency domain signal-to-noise ratio, more accurate echo signal broadband time delay vector is calculated, the reconstruct source signal of fidelity enhancing is obtained.
Description
Technical field
The invention belongs to signal processing fields more particularly to a kind of broadband time delay for high-fidelity ARRAY PROCESSING to estimate
Meter method.
Background technique
Acoustic Object Classification and Identification is one of important topic of Underwater acoustic signal processing.Passive sonar target classification is known
Not, it is important to the feature extraction of target radiated noise.In order to improve the performance of feature extraction, Wave beam forming processing is from sensor
Source signal is reconstructed in array signal, the effective means and common method of spatially antinoise and separation multi-target jamming.But
In Underwater Acoustic Environment, due to the error of element position, the hypothesis of plane wave and evaluated error of azimuth of target etc., may all it make
There is delay mismatch when exporting at conventional beamformer method tracking beam, to cause signal special in power spectrum and demodulation spectra etc.
The distortion of aspect is levied, to influence the effect of target classification identification.So this special sonar mission is appointed for target identification
Business, it is the research topic that must pay attention to that high-fidelity, which reconstructs source signal and the broadband array processing method of feature,.
Domestic and foreign scholars propose many high-fidelity broadband array processing methods at present, mainly have: (1) array calibration: passing through
Evaluated error parameter realizes the correction to formation, such method first models array error, array error is corrected
The problem of being converted into parameter Estimation.Such array calibration method can be generally divided into active correction class and self-correcting class.For having
Source bearing calibration, this method have the requirement of higher refined orientation information to auxiliary source, so when auxiliary source
When azimuth information has deviation, this kind of correcting algorithm can bring biggish deviation.Self-Tuning Algorithm is due to sensor position uncertainties and side
The array structure of coupling and certain morbid state between the parameter of position, the unique identification of parameter Estimation often not can guarantee, even more important
Be that parametric joint estimates that corresponding higher-dimension, multimode nonlinear optimal problem bring huge operand, the global of estimation is received
Holding back property often not can guarantee.(2) spherical wave focused beamforming: is carried out to delay inequality according to different direction in measurement zone and distance
Compensation obtains the distribution map of target sound source in measurement zone, to realize the high-precision Passive Positioning of near-field target, according to the mesh of estimation
Using conventional spherical wave Wave beam forming high-fidelity reconstruct source signal, the technical principle is simple, is easily achieved, energy for subject distance and orientation
It is enough effectively to solve the problems, such as under the conditions of near field sources delay mismatch between plane wave approximation bring array element, but operand greatly and not can solve
Sensor position uncertainties bring delay mismatch (3) ESPRIT (estimating signal parameter by rotational invariance) super-resolution direction finding: benefit
Parameter Estimation is carried out with the translation invariance of array, signal parameter directly can be obtained by characteristic value, this method is without carrying out
The great spectrum peak search of calculation amount is not necessarily to storage array prevalence matrix, but requires accurately known (4) the constant beam-width wave beam shape of formation
At: the weighing vector of selection Wave beam forming in each frequency component in broadband signal bandwidth, so that basic matrix system is to different frequency
Input signal space filtering having the same response, which can eliminate tracking beam caused by azimuth of target evaluated error
The non-linear distortion of output signal, but operand it is big and not can solve due to the error of element position, the hypothesis of plane wave and
(5) zero pilot beam of delay mismatch caused by evaluated error of azimuth of target etc. is formed: adjustment weight makes battle array response in signal
Direction is 1, is 0 in interference radiating way, to achieve the purpose that inhibit interference, it is actually to form beam pattern in interference radiating way
Zero point, which plays an important role to directional jamming is eliminated, but needs to signal and jamming bearing and formation information accurately
Know.
Summary of the invention
Goal of the invention: aiming at the problems existing in the prior art, the invention discloses one kind to be used for high-fidelity ARRAY PROCESSING
Broadband delay time estimation method, this method is based on time delay between the power spectrum line spectrum array element having the same of echo signal different frequency
This characteristic is added phase difference between the array element estimated based on different frequency line spectrum according to the frequency of line spectrum and frequency domain signal-to-noise ratio
Power fusion, calculates more accurate echo signal broadband time delay vector, obtains the reconstruct source signal of fidelity enhancing.
Technical solution: in order to achieve the above-mentioned object of the invention, the invention adopts the following technical scheme: a kind of be used for high-fidelity battle array
The broadband delay time estimation method of column processing, includes the following steps:
(1) based on echo signal far field plane wave incidence and undistorted two of formation it is assumed that calculating basic matrix acquisition data not
Pre- with orientation seeks beam energy figure at beam data, estimates target bearing
(2) based on echo signal far field plane wave incidence and undistorted two of formation it is assumed that being pressed from basic matrix acquisition data
The target bearing of estimationSource signal is reconstructed in advance, detects pre- reconstruct source signal s1[n], n=1, the power spectrum of 2 ..., N
Line spectrum position fk, k=1,2 ..., K calculate the accurate measurement value of line spectral frequenciesWherein K is the frequency detected
Power spectrum line spectrum radical less than V/ (2d), V are the spread speed of sound wave in water, and d is even linear array array element spacing；
(3) reference array element H is calculated1The power spectrum of acquisition data existsThe frequency domain Signal to Noise Ratio (SNR) of line spectrum at positionk, k=1,
2 ..., K；
(4) frequency spectrum for calculating basic matrix acquisition data existsPhase difference between the array element of line spectrum at position Wherein, M is basic matrix array element number；
(5) according to line spectral frequencies and frequency domain signal-to-noise ratio to phase difference between the array element estimated based on different frequency line spectrum
It is weighted fusion, estimates echo signal broadband time delay vector
(6) according to the echo signal broadband time delay vector of estimationSource signal is carried out from basic matrix acquisition data secondary heavy
Structure obtains the reconstruct source signal s of fidelity enhancing2[n], n=1,2 ..., N.
Preferably, in step (1), the beam energy figure of basic matrix acquisition data is calculated with the following method, estimates target side
Position, specifically comprises the following steps:
(1-1) calculates pre- at beam angle by waiting cosine to divide:
Wherein, θiFlat square is established using first array element of basic matrix as coordinate origin at the angle of wave beam in advance for i-th
Coordinate system, the coordinate of m-th of array element of basic matrix are [- (m-1) d, 0], m=1,2 ..., M, θiAngle is defined as origin and target is sat
The angle of target line and x-axis, θi∈ [0, π],It indicates to calculateAnticosine, I be in advance at
Wave beam number, d are even linear array array element spacing；
(1-2) is based on echo signal far field plane wave incidence and undistorted two of formation it is assumed that calculating even linear array is adjacent
Array element acquires signal in advance at beam angle θiUnder delay inequality τi:
Wherein, cos (θi) indicate to calculate θiCosine value, V be the spread speed of sound wave in water；
(1-3) presses following beamforming algorithm；With first array element H of basic matrix battle array head1For reference array element, calculate in advance at wave beam
Angle, θiIt is corresponding pre- at beam data
Wherein, M is basic matrix array element number, xm[n] indicates m-th of array element, n-th of sampled data, and N is in advance at beam data
Discrete points, fsSample frequency when data is acquired for basic matrix；
(1-4) calculates basic matrix acquisition data in advance at beam angle θiUnder beam energy E [i]:
(1-5) estimates target bearing
It searches for pre- at beam index i corresponding to E [i] maximum valueMax
Wherein, argmax1≤i≤I{ E [i] } indicates to search for pre- Cheng Bo corresponding to E [i] maximum value within the scope of 1≤i≤I
Beam indexes iMax.If 1 < iMax< I, willIt is denoted as Z respectively-1、Z0、Z1, then base
In the target bearing of parabola three point interpolation estimationAre as follows:
Preferably, in step (2), detection reconstructs the power spectrum line spectrum position f of source signal in advance with the following methodk, meter
Calculate the accurate measurement value of line spectral frequenciesSpecifically comprise the following steps:
(2-1) presses following beamforming algorithm, with first array element H of basic matrix battle array head1For reference array element, number is acquired from basic matrix
By the target bearing of estimation inSource signal is reconstructed in advance, obtains pre- reconstruct source signal s1[n]:
(2-2) calculates pre- reconstruct source signal s1The discrete spectrum of [n]:
Source signal s is then reconstructed in advance1The power spectrum of [n] are as follows:
Wherein, l S1[l] and P1The discrete frequency of [l] indexes, and j indicates imaginary unit, i.e.,| | represent modulus
It is worth operation；
(2-3) is by power spectrum P1[l] obtains power spectrum P by Q rank median filter1The trend term C of [l]1[l], l=0,
1 ..., N/2, then power spectrum P1[l] and trend term C1The difference of [l] composes D1[l] are as follows:
D1[l]=P1[l]-C1[l], l=0,1 ..., N/2
D1The standardization difference of [l] is composed are as follows:
Wherein, std (D1[l]) it is to seek D1The standard deviation of [l]；
(2-4) sets amplitude threshold G and extracts power spectrum line spectrum, if D2(l) meet following condition and be then judged to power spectrum line spectrum:
Assuming that extracting the power spectrum line spectrum that K root frequency is less than V/2d altogether, kth root line spectrum is in D2Index in [l] is Jk,
The then frequency of kth root line spectrum are as follows:
fk=JkΔ f, k=1,2 ..., K
Wherein, Δ f is the frequency resolution for the discrete Fourier transform that length is N, Δ f=fs/N；
(2-5) takes S1[l] indexes (J in discrete frequencyk- 1), Jk, (Jk+ 1) modulus value at, is denoted as A respectivelyk1, Ak2, Ak3,
That is Ak1=| S1[Jk- 1] |, Ak2=| S1[Jk] |, Ak3=| S1[Jk+ 1] |, utilize the modulus value Ak1, Ak2, Ak3Calculate kth root function
The relative deviation δ of rate spectral line spectral frequency interpolationk, it may be assumed that
The then frequency accurate measurement value of kth root power spectrum line spectrumAre as follows:
Preferably, in step (3), reference array element H is calculated with the following method1Acquire data power spectrumLine at position
The frequency domain Signal to Noise Ratio (SNR) of spectrumk, specifically comprise the following steps:
(3-1) calculates reference array element H1Acquire data x1The power spectrum of [n]:
Wherein, l isDiscrete frequency index；
(3-2) estimates reference array element H1Acquire data x1The power spectrum of [n] existsThe performance number of noise at position
Wherein, WL is the length of noise power-value estimator, takes even number, and WL≤2 (Jk-1)；
(3-3) calculates reference array element H1Acquire data x1The power spectrum of [n] existsThe performance number of signal plus noise at position
(3-4) reference array element H1Acquire data x1The power spectrum of [n] existsThe frequency domain Signal to Noise Ratio (SNR) of line spectrum at positionkAre as follows:
Preferably, in step (4), the frequency spectrum for calculating basic matrix acquisition data with the following method existsArray element at position
Between phase differenceSpecifically comprise the following steps:
(4-1) calculates m-th of array element HmThe frequency spectrum of acquisition data existsValue at position:
The frequency spectrum that (4-2) calculates each array element acquisition data existsPhase vectors at position
Wherein,To seekThe operation of phase；
The frequency spectrum that (4-3) calculates adjacent array element acquisition data existsThe difference of phase at position
Wherein,For vectorM-th of element, m=1,2 ..., M；
Preferably, in step (5), with the following method between the echo signal array element estimated based on different frequency line spectrum
Phase difference is weighted fusion, estimates echo signal broadband time delay vectorSpecifically comprise the following steps:
(5-1) is according to line spectral frequenciesWith frequency domain Signal to Noise Ratio (SNR)k, calculate phase difference weighting system between the array element of kth root line spectrum
Number wk:
(5-2) is according to phase difference weighting coefficient w between array elementkBetween phase difference array elementIt is weighted fusion, estimates mesh
Mark signal broadband time delay vector
Wherein
Preferably, in step (6), secondary reconstruct is carried out to source signal with the following method: being believed according to the target of estimation
Number broadband time delay vectorData x is acquired from basic matrixmSecondary reconstruct is carried out to source signal in [n], obtains the reconstruct source of fidelity enhancing
Signal s2[n]:
The utility model has the advantages that compared with the existing methods, technical solution of the present invention has following advantageous effects:
(1) method disclosed by the invention estimates echo signal broadband based on phase difference between the array element of the strong line spectrum of echo signal
Time delay vector realizes formation distortion, data-driven, adaptive fidelity enhancing source signal reconstruct under the conditions of near field sources；
(2) method of the invention based on time delay between the power spectrum line spectrum array element having the same of echo signal different frequency this
One characteristic, according to the frequency of power spectrum line spectrum and frequency domain signal-to-noise ratio to phase difference between the array element estimated based on different frequency line spectrum into
Row Weighted Fusion comprehensively utilizes phase information between the array element of each line spectrum of target, calculates more accurate echo signal broadband time delay
Vector obtains the reconstruct source signal of fidelity enhancing.
Detailed description of the invention
Fig. 1 is the flow chart of the method for the present invention；
Fig. 2 is the array element coordinate position comparison diagram of distortion battle array and undistorted battle array in embodiment；
Fig. 3 is target and each element position schematic diagram in embodiment；
Fig. 4 be in embodiment known to formation and target range when beam energy figure be based on echo signal far field plane wave
Incident and the undistorted hypothesis of formation beam energy figure；
Fig. 5 is the power spectrum that reference array element acquires data in embodiment；
Fig. 6 is digital time delay comparison diagram between the echo signal array element calculated in different ways in embodiment；
Fig. 7 is the power spectrum line spectrum feature comparison diagram that source signal is reconstructed in embodiment；
Fig. 8 is the demodulation spectra line spectrum feature comparison diagram that source signal is reconstructed in embodiment.
Specific embodiment
Further description of the technical solution of the present invention with reference to the accompanying drawings and examples.
As shown in Figure 1, a kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING, includes the following steps:
(1) based on echo signal far field plane wave incidence and undistorted two of formation it is assumed that calculating basic matrix acquisition data not
Pre- with orientation seeks beam energy figure at beam data, estimates target bearing
(2) based on echo signal far field plane wave incidence and undistorted two of formation it is assumed that being pressed from basic matrix acquisition data
The target bearing of estimationSource signal is reconstructed in advance, detects pre- reconstruct source signal s1[n], n=1, the power spectrum of 2 ..., N
Line spectrum position fk, k=1,2 ..., K calculate the accurate measurement value of line spectral frequenciesWherein K is the frequency detected
Power spectrum line spectrum radical less than V/ (2d), V are the spread speed of sound wave in water, and d is even linear array array element spacing；
(3) reference array element H is calculated1The power spectrum of acquisition data existsThe frequency domain Signal to Noise Ratio (SNR) of line spectrum at positionk, k=1,
2 ..., K；
(4) frequency spectrum for calculating basic matrix acquisition data existsPhase difference between the array element of line spectrum at position Wherein, M is basic matrix array element number；
(5) according to line spectral frequencies and frequency domain signal-to-noise ratio to phase difference between the array element estimated based on different frequency line spectrum
It is weighted fusion, estimates echo signal broadband time delay vector
(6) according to the echo signal broadband time delay vector of estimationSource signal is carried out from basic matrix acquisition data secondary heavy
Structure obtains the reconstruct source signal s of fidelity enhancing2[n], n=1,2 ..., N.
In step (1), the beam energy figure of basic matrix acquisition data is calculated with the following method, estimates target bearing, specifically
Include the following steps:
(1-1) calculates pre- at beam angle by waiting cosine to divide:
Wherein, θiFlat square is established using first array element of basic matrix as coordinate origin at the angle of wave beam in advance for i-th
Coordinate system, the coordinate of m-th of array element of basic matrix are [- (m-1) d, 0], and m=1,2 ..., M, target angle are defined as origin and target
The line of coordinate and the angle of x-axis, θi∈ [0, π],It indicates to calculateAnticosine, I is pre-
At wave beam number, d is even linear array array element spacing；
(1-2) is based on echo signal far field plane wave incidence and undistorted two of formation it is assumed that calculating even linear array is adjacent
Array element acquires signal in advance at beam angle θiUnder delay inequality τi:
Wherein, cos (θi) indicate to calculate θiCosine value, V be the spread speed of sound wave in water；
(1-3) presses following beamforming algorithm；With first array element H of basic matrix battle array head1For reference array element, calculate in advance at wave beam
Angle, θiIt is corresponding pre- at beam data
Wherein, M is basic matrix array element number, xm[n] indicates m-th of array element, n-th of sampled data, and N is in advance at beam data
Discrete points, fsSample frequency when data is acquired for basic matrix；
(1-4) calculates basic matrix acquisition data in advance at beam angle θiUnder beam energy E [i]:
(1-5) estimates target bearing
It searches for pre- at beam index i corresponding to E [i] maximum valueMax
Wherein, argmax1≤i≤I{ E [i] } indicates to search for pre- Cheng Bo corresponding to E [i] maximum value within the scope of 1≤i≤I
Beam indexes iMax.If 1 < iMax< I, willIt is denoted as Z respectively-1、Z0、Z1, then base
In the target bearing of parabola three point interpolation estimationAre as follows:
In step (2), detection reconstructs the power spectrum line spectrum position f of source signal in advance with the following methodk, calculate line spectrum frequency
The accurate measurement value of rateSpecifically comprise the following steps:
(2-1) presses following beamforming algorithm, with first array element H of basic matrix battle array head1For reference array element, number is acquired from basic matrix
By the target bearing of estimation inSource signal is reconstructed in advance, obtains pre- reconstruct source signal s1[n]:
(2-2) calculates pre- reconstruct source signal s1The discrete spectrum of [n]:
Source signal s is then reconstructed in advance1The power spectrum of [n] are as follows:
Wherein, l S1[l] and P1The discrete frequency of [l] indexes, and j indicates imaginary unit, i.e.,| | represent modulus
It is worth operation；
(2-3) is by power spectrum P1[l] obtains power spectrum P by Q rank median filter1The trend term C of [l]1[l], l=0,
1 ..., N/2, then power spectrum P1[l] and trend term C1The difference of [l] composes D1[l] are as follows:
D1[l]=P1[l]-C1[l], l=0,1 ..., N/2
D1The standardization difference of [l] is composed are as follows:
Wherein, std (D1[l]) it is to seek D1The standard deviation of [l]；
(2-4) sets amplitude threshold G and extracts power spectrum line spectrum, if D2(l) meet following condition and be then judged to power spectrum line spectrum:
Assuming that extracting the power spectrum line spectrum that K root frequency is less than V/2d altogether, kth root line spectrum is in D2Index in [l] is Jk,
The then frequency of kth root line spectrum are as follows:
fk=JkΔ f, k=1,2 ..., K
Wherein, Δ f is the frequency resolution for the discrete Fourier transform that length is N, Δ f=fs/N；
(2-5) takes S1[l] indexes (J in discrete frequencyk- 1), Jk, (Jk+ 1) modulus value at, is denoted as A respectivelyk1, Ak2, Ak3,
That is Ak1=| S1[Jk- 1] |, Ak2=| S1[Jk] |, Ak3=| S1[Jk+ 1] |, utilize the modulus value Ak1, Ak2, Ak3Calculate kth root function
The relative deviation δ of rate spectral line spectral frequency interpolationk, it may be assumed that
The then frequency accurate measurement value of kth root power spectrum line spectrumAre as follows:
In step (3), reference array element H is calculated with the following method1Acquire data power spectrumThe frequency domain of line spectrum at position
Signal to Noise Ratio (SNR)k, specifically comprise the following steps:
(3-1) calculates reference array element H1Acquire data x1The power spectrum of [n]:
Wherein, l isDiscrete frequency index；
(3-2) estimates reference array element H1Acquire data x1The power spectrum of [n] existsThe performance number of noise at position
Wherein, WL is the length of noise power-value estimator, takes even number, and WL≤2 (Jk-1)；
(3-3) calculates reference array element H1Acquire data x1The power spectrum of [n] existsThe performance number of signal plus noise at position
(3-4) reference array element H1Acquire data x1The power spectrum of [n] existsThe frequency domain Signal to Noise Ratio (SNR) of line spectrum at positionkAre as follows:
Preferably, in step (4), the frequency spectrum for calculating basic matrix acquisition data with the following method existsArray element at position
Between phase differenceSpecifically comprise the following steps:
(4-1) calculates m-th of array element HmThe frequency spectrum of acquisition data existsValue at position:
The frequency spectrum that (4-2) calculates each array element acquisition data existsPhase vectors at position
Wherein,To seekThe operation of phase；
The frequency spectrum that (4-3) calculates adjacent array element acquisition data existsThe difference of phase at position
Wherein,For vectorM-th of element, m=1,2 ..., M；
In step (5), with the following method between phase difference the echo signal array element estimated based on different frequency line spectrum
It is weighted fusion, estimates echo signal broadband time delay vectorSpecifically comprise the following steps:
(5-1) is according to line spectral frequenciesWith frequency domain Signal to Noise Ratio (SNR)k, calculate phase difference weighting system between the array element of kth root line spectrum
Number wk:
(5-2) is according to phase difference weighting coefficient w between array elementkBetween phase difference array elementIt is weighted fusion, estimates mesh
Mark signal broadband time delay vector
Wherein
In step (6), secondary reconstruct is carried out to source signal with the following method: when according to the echo signal broadband of estimation
Prolong vectorData x is acquired from basic matrixmSecondary reconstruct is carried out to source signal in [n], obtains the reconstruct source signal s of fidelity enhancing2
[n]:
Embodiment:
Emulation signal parameter is respectively set are as follows: even linear array array element number M=32, array element spacing d=0.375m, distort battle array
It is as shown in Figure 2 with undistorted a burst of first coordinate position comparison.Target and undistorted a burst of first coordinate position signal as shown in figure 3,
Using first first array element of battle array as coordinate origin, angle theta=120 ° of target and the direction battle array head, target is at a distance from reference array element
R=3000m.Sample frequency fs=10kHz, the propagation speed V=1500m/s of sound wave in water.Single array element receives data SNR
For -9.62dB, comprising 100.08Hz, 300.1Hz, 500.2Hz, 800.3Hz, 1000.41Hz, 1200.7Hz, 1500.5Hz,
1801.5Hz, 1999.5Hz totally 9 line spectrums, when a length of 8s of time analysis window, corresponding frequency domain signal-to-noise ratio be respectively -3dB,
12dB,0dB,16dB,18dB,-3dB,-3dB,-3dB,10dB；2Hz, 4Hz, 6Hz, 8Hz totally 4 modulation spectrum line spectrums are set, are adjusted
Depth processed is respectively 0.04,0.04,0.04,0.06.
In (1) step, setting in advance at wave beam number I=65, based on echo signal far field plane wave incidence and formation without
Distortion two it is assumed that the beam energy figure sought as shown in figure 4, the target bearing of estimation is
According to (2) step, the order that median filter is arranged is Q=51, and standardization difference composes extraction of line spectrum amplitude threshold G
=20.The power spectrum line spectrum radical detected is K=4, and the frequency of 4 power spectrum line spectrums is respectively
The frequency accurate measurement value of 4 power spectrum line spectrums based on Interpolate estimation is respectively
According to (3) step, reference array element acquires the power spectrum of data as shown in figure 5, when calculating line spectrum frequency domain signal-to-noise ratio, joins
With average frequency points WL=80, array element H is examined1The frequency domain signal-to-noise ratio of power spectrum respective frequencies line spectrum for acquiring data is respectively
According to (5) step, phase difference weighting coefficient is respectively between the array element of 4 power spectrum line spectral frequencies
Fusion is weighted between phase difference the array element of 4 line spectrums according to weighting coefficient, estimates echo signal broadband time delay
VectorDigital time delay is as shown in Figure 6 between the echo signal array element calculated in different ways, it can be seen that relative to based on target
Digital time delay between signal far field plane wave incidence and the undistorted array element for assuming estimation of formation, is based on the multi-thread spectrum array element of echo signal
Between phase difference partial combination assessing array element between digital time delay and echo signal array element digital time delay true value have better one
Cause property.
According to (6) step, according to the echo signal broadband time delay vector of estimationTo source signal from basic matrix acquisition data
Secondary reconstruct is carried out, the reconstruct source signal of fidelity enhancing is obtained.Fig. 7 is pre- reconstruct source signal power spectrum and secondary reconstruct source signal
Power spectrum line spectrum feature comparison diagram, it can be seen that compared to pre- reconstruct source signal, the power spectrum line spectrum of secondary reconstruct source signal is special
Sign fidelity is remarkably reinforced.For the power spectrum line spectrum that can be detected from pre- reconstruct source signal, examined from secondary reconstruct source signal
The power spectrum line spectrum measured has higher signal-to-noise ratio, and can detect from secondary reconstruct source signal can not be from pre- reconstruct source signal
The power spectrum line spectrum detected.Fig. 8 is that pre- reconstruct source signal demodulation spectra and secondary reconstruct source signal demodulation spectra line spectrum feature compare
Figure, it can be seen that compared to pre- reconstruct source signal, the demodulation spectra structure of secondary reconstruct source signal is apparent, line spectrum feature noise
Than higher.
Claims (7)
1. a kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING, which is characterized in that this method comprises the following steps:
(1) based on echo signal far field plane wave incidence and undistorted two of formation it is assumed that calculating basic matrix acquisition data not Tongfang
The pre- of position seeks beam energy figure at beam data, estimates target bearing
(2) based on echo signal far field plane wave incidence and undistorted two of formation it is assumed that by estimation from basic matrix acquisition data
Target bearingSource signal is reconstructed in advance, detects pre- reconstruct source signal s1The power spectrum line spectrum position f of [n]k, calculate line
The accurate measurement value of spectral frequencyWherein, k=1,2 ..., K, K are that the frequency detected is less than the power spectrum line spectrum radical of V/ (2d),
V is the spread speed of sound wave in water, and d is even linear array array element spacing, and N is in advance the n=1 at beam data discrete points,
2 ..., N；
(3) reference array element H is calculated1The power spectrum of acquisition data existsThe frequency domain Signal to Noise Ratio (SNR) of line spectrum at positionk；
(4) frequency spectrum for calculating basic matrix acquisition data existsPhase difference between the array element of line spectrum at position
Wherein, M is basic matrix array element number；
(5) according to line spectral frequencies and frequency domain signal-to-noise ratio to phase difference between the array element estimated based on different frequency line spectrumAdded
Power fusion, estimates echo signal broadband time delay vector
(6) according to the echo signal broadband time delay vector of estimationSecondary reconstruct is carried out to source signal from basic matrix acquisition data, is obtained
Go bail for the reconstruct source signal s really enhanced2[n]。
2. a kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING according to claim 1, which is characterized in that
In step (1), the beam energy figure of basic matrix acquisition data is calculated with the following method, estimates target bearingIt specifically includes as follows
Step:
(1-1) calculates pre- at beam angle by waiting cosine to divide:
Wherein, θiPlane rectangular coordinates is established using first array element of basic matrix as coordinate origin at the angle of wave beam in advance for i-th
System, the coordinate of m-th of array element of basic matrix are [- (m-1) d, 0], m=1,2 ..., M, θiAngle is defined as origin and coordinates of targets
The angle of line and x-axis, θi∈ [0, π],It indicates to calculateAnticosine, I be pre- Cheng Bo
Beam number, d are even linear array array element spacing；
(1-2) is based on echo signal far field plane wave incidence and undistorted two of formation it is assumed that calculating the adjacent array element of even linear array
Signal is acquired in advance at beam angle θiUnder delay inequality τi:
Wherein, cos (θi) indicate to calculate θiCosine value, V be the spread speed of sound wave in water；
(1-3) presses following beamforming algorithm, with first array element H of basic matrix battle array head1For reference array element, calculate in advance at beam angle
θiIt is corresponding pre- at beam data
Wherein, M is basic matrix array element number, xm[n] indicates m-th of array element, n-th of sampled data, and N is in advance at beam data discrete point
Number, fsSample frequency when data is acquired for basic matrix；
(1-4) calculates basic matrix acquisition data in advance at beam angle θiUnder beam energy E [i]:
(1-5) estimates target bearing
It searches for pre- at beam index i corresponding to E [i] maximum valueMax
Wherein, argmax1≤i≤I{ E [i] } indicates to search within the scope of 1≤i≤I pre- at wave beam rope corresponding to E [i] maximum value
Draw iMaxIf 1 < iMax< I, willIt is denoted as Z respectively-1、Z0、Z1, then it is based on
The target bearing of parabola three point interpolation estimationAre as follows:
3. a kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING according to claim 2, which is characterized in that
In step (2), detection reconstructs the power spectrum line spectrum position f of source signal in advance with the following methodk, calculate the accurate measurement of line spectral frequencies
ValueSpecifically comprise the following steps:
(2-1) presses following beamforming algorithm, with first array element H of basic matrix battle array head1For reference array element, acquired in data from basic matrix
By the target bearing of estimationSource signal is reconstructed in advance, obtains pre- reconstruct source signal s1[n]:
(2-2) calculates pre- reconstruct source signal s1The discrete spectrum of [n]:
Source signal s is then reconstructed in advance1The power spectrum of [n] are as follows:
Wherein, l S1[l] and P1The discrete frequency of [l] indexes, and j indicates imaginary unit, i.e.,| | represent modulus value fortune
It calculates；
(2-3) is by power spectrum P1[l] obtains power spectrum P by Q rank median filter1The trend term C of [l]1[l], l=0,
1 ..., N/2, then power spectrum P1[l] and trend term C1The difference of [l] composes D1[l] are as follows:
D1[l]=P1[l]-C1[l], l=0,1 ..., N/2
D1The standardization difference of [l] is composed are as follows:
Wherein, std (D1[l]) it is to seek D1The standard deviation of [l]；
(2-4) sets amplitude threshold G and extracts power spectrum line spectrum, if D2(l) meet following condition and be then judged to power spectrum line spectrum:
Assuming that extracting the power spectrum line spectrum that K root frequency is less than V/2d altogether, kth root line spectrum is in D2Index in [l] is Jk, then kth
The frequency of root line spectrum are as follows:
fk=JkΔ f, k=1,2 ..., K
Wherein, Δ f is the frequency resolution for the discrete Fourier transform that length is N, Δ f=fs/N；
(2-5) takes S1[l] indexes (J in discrete frequencyk- 1), Jk, (Jk+ 1) modulus value at, is denoted as A respectivelyk1, Ak2, Ak3, i.e. Ak1
=| S1[Jk- 1] |, Ak2=| S1[Jk] |, Ak3=| S1[Jk+ 1] |, utilize the modulus value Ak1, Ak2, Ak3Calculate kth root power spectrum
The relative deviation δ of line spectral frequencies interpolationk, it may be assumed that
The then frequency accurate measurement value of kth root power spectrum line spectrumAre as follows:
4. a kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING according to claim 3, which is characterized in that
In step (3), reference array element H is calculated with the following method1Acquire data power spectrumThe frequency domain signal-to-noise ratio of line spectrum at position
SNRk, specifically comprise the following steps:
(3-1) calculates reference array element H1Acquire data x1The power spectrum of [n]:
Wherein, l isDiscrete frequency index；
(3-2) estimates reference array element H1Acquire data x1The power spectrum of [n] existsThe performance number of noise at position
Wherein, WL is the length of noise power-value estimator, takes even number, and WL≤2 (Jk-1)；
(3-3) calculates reference array element H1Acquire data x1The power spectrum of [n] existsThe performance number of signal plus noise at position
(3-4) reference array element H1Acquire data x1The power spectrum of [n] existsThe frequency domain Signal to Noise Ratio (SNR) of line spectrum at positionkAre as follows:
5. a kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING according to claim 4, which is characterized in that
In step (4), the frequency spectrum for calculating basic matrix acquisition data with the following method existsPhase difference between array element at position
Specifically comprise the following steps:
(4-1) calculates m-th of array element HmThe frequency spectrum of acquisition data existsValue at position:
The frequency spectrum that (4-2) calculates each array element acquisition data existsPhase vectors at position
Wherein,To seekThe operation of phase；
The frequency spectrum that (4-3) calculates adjacent array element acquisition data existsThe difference of phase at position
Wherein,For vectorM-th of element, m=1,2 ..., M.
6. a kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING according to claim 5, which is characterized in that
In step (5), phase difference the echo signal array element estimated based on different frequency line spectrum is weighted with the following method
Echo signal broadband time delay vector is estimated in fusionSpecifically comprise the following steps:
(5-1) is according to line spectral frequenciesWith frequency domain Signal to Noise Ratio (SNR)k, calculate phase difference weighting coefficient w between the array element of kth root line spectrumk:
(5-2) is according to phase difference weighting coefficient w between array elementkBetween phase difference array elementIt is weighted fusion, estimation target letter
Number broadband time delay vector
Wherein,
7. a kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING according to claim 6, which is characterized in that
In step (6), secondary reconstruct is carried out to source signal with the following method: according to the echo signal broadband time delay vector of estimation
Data x is acquired from basic matrixmSecondary reconstruct is carried out to source signal in [n], obtains the reconstruct source signal s of fidelity enhancing2[n]:
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910001466.6A CN109799495B (en) | 2019-01-02 | 2019-01-02 | Broadband time delay estimation method for high fidelity array processing |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910001466.6A CN109799495B (en) | 2019-01-02 | 2019-01-02 | Broadband time delay estimation method for high fidelity array processing |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109799495A true CN109799495A (en) | 2019-05-24 |
| CN109799495B CN109799495B (en) | 2022-06-17 |
Family
ID=66558408
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910001466.6A Active CN109799495B (en) | 2019-01-02 | 2019-01-02 | Broadband time delay estimation method for high fidelity array processing |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109799495B (en) |
Cited By (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110208808A (en) * | 2019-05-28 | 2019-09-06 | 哈尔滨工程大学 | A kind of passive sonar noncooperative target line spectrum information fusion method |
| CN111487478A (en) * | 2020-03-27 | 2020-08-04 | 杭州电子科技大学 | Angle-dependent complex array error calibration method based on deep neural network |
| CN111929666A (en) * | 2020-09-09 | 2020-11-13 | 东南大学 | Weak underwater sound target line spectrum autonomous extraction method based on sequential environment learning |
| CN112379357A (en) * | 2020-11-16 | 2021-02-19 | 东南大学 | Signal discrimination method based on pulse signal pre-estimation time parameter secondary processing |
| CN112816938A (en) * | 2020-12-25 | 2021-05-18 | 中国船舶重工集团有限公司第七一0研究所 | Active sonar pulse signal direction finding method |
| CN113109760A (en) * | 2021-04-13 | 2021-07-13 | 东南大学 | Multi-line spectrum combined DOA estimation and clustering method and system based on group sparsity |
| CN113325401A (en) * | 2021-07-06 | 2021-08-31 | 东南大学 | Distortion towed linear array signal reconstruction method based on line spectrum phase difference ambiguity resolution |
| CN113625222A (en) * | 2021-08-27 | 2021-11-09 | 东南大学 | Broadband signal high-resolution direction estimation method based on time delay difference norm weighting |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102621527A (en) * | 2012-03-20 | 2012-08-01 | 哈尔滨工程大学 | Broad band coherent source azimuth estimating method based on data reconstruction |
| CN104101871A (en) * | 2013-04-15 | 2014-10-15 | 中国科学院声学研究所 | Narrowband interference suppression method and narrowband interference suppression system used for passive synthetic aperture |
| US20150073276A1 (en) * | 2005-07-22 | 2015-03-12 | Zonare Medical Systems, Inc. | Aberration correction using channel data in ultrasound imaging system |
| CN105022050A (en) * | 2014-04-16 | 2015-11-04 | 中国科学院声学研究所 | Underwater-acoustic-channel discrete noise source suppression method of multi-sensor array |
| CN107132503A (en) * | 2017-03-23 | 2017-09-05 | 哈尔滨工程大学 | Acoustic vector circle battle array broadband coherent source direction estimation method based on vector singular value decomposition |
| CN108761383A (en) * | 2018-04-13 | 2018-11-06 | 中国人民解放军陆军工程大学 | A kind of time delay based on two-dimensional matrix beam and angle combined estimation method |
-
2019
- 2019-01-02 CN CN201910001466.6A patent/CN109799495B/en active Active
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20150073276A1 (en) * | 2005-07-22 | 2015-03-12 | Zonare Medical Systems, Inc. | Aberration correction using channel data in ultrasound imaging system |
| CN102621527A (en) * | 2012-03-20 | 2012-08-01 | 哈尔滨工程大学 | Broad band coherent source azimuth estimating method based on data reconstruction |
| CN104101871A (en) * | 2013-04-15 | 2014-10-15 | 中国科学院声学研究所 | Narrowband interference suppression method and narrowband interference suppression system used for passive synthetic aperture |
| CN105022050A (en) * | 2014-04-16 | 2015-11-04 | 中国科学院声学研究所 | Underwater-acoustic-channel discrete noise source suppression method of multi-sensor array |
| CN107132503A (en) * | 2017-03-23 | 2017-09-05 | 哈尔滨工程大学 | Acoustic vector circle battle array broadband coherent source direction estimation method based on vector singular value decomposition |
| CN108761383A (en) * | 2018-04-13 | 2018-11-06 | 中国人民解放军陆军工程大学 | A kind of time delay based on two-dimensional matrix beam and angle combined estimation method |
Non-Patent Citations (2)
| Title |
|---|
| 方世良等: ""基于频域分解的水声信号宽带盲波束形成"", 《声学技术》 * |
| 范展 等: ""基于干扰协方差矩阵重构的恒定束宽鲁棒自适应波束形成"", 《声学学报》 * |
Cited By (13)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110208808B (en) * | 2019-05-28 | 2022-07-15 | 哈尔滨工程大学 | Passive sonar non-cooperative target line spectrum information fusion method |
| CN110208808A (en) * | 2019-05-28 | 2019-09-06 | 哈尔滨工程大学 | A kind of passive sonar noncooperative target line spectrum information fusion method |
| CN111487478A (en) * | 2020-03-27 | 2020-08-04 | 杭州电子科技大学 | Angle-dependent complex array error calibration method based on deep neural network |
| CN111929666A (en) * | 2020-09-09 | 2020-11-13 | 东南大学 | Weak underwater sound target line spectrum autonomous extraction method based on sequential environment learning |
| CN112379357A (en) * | 2020-11-16 | 2021-02-19 | 东南大学 | Signal discrimination method based on pulse signal pre-estimation time parameter secondary processing |
| CN112379357B (en) * | 2020-11-16 | 2023-10-27 | 东南大学 | Signal discrimination method based on secondary processing of pulse signal estimated time parameter |
| CN112816938B (en) * | 2020-12-25 | 2022-08-09 | 中国船舶重工集团有限公司第七一0研究所 | Active sonar pulse signal direction finding method |
| CN112816938A (en) * | 2020-12-25 | 2021-05-18 | 中国船舶重工集团有限公司第七一0研究所 | Active sonar pulse signal direction finding method |
| CN113109760A (en) * | 2021-04-13 | 2021-07-13 | 东南大学 | Multi-line spectrum combined DOA estimation and clustering method and system based on group sparsity |
| CN113325401A (en) * | 2021-07-06 | 2021-08-31 | 东南大学 | Distortion towed linear array signal reconstruction method based on line spectrum phase difference ambiguity resolution |
| CN113325401B (en) * | 2021-07-06 | 2024-03-19 | 东南大学 | Distortion towing linear array signal reconstruction method based on line spectrum phase difference deblurring |
| CN113625222A (en) * | 2021-08-27 | 2021-11-09 | 东南大学 | Broadband signal high-resolution direction estimation method based on time delay difference norm weighting |
| CN113625222B (en) * | 2021-08-27 | 2022-10-28 | 东南大学 | Broadband signal high-resolution direction estimation method based on time delay difference norm weighting |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109799495B (en) | 2022-06-17 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109799495A (en) | A kind of broadband delay time estimation method for high-fidelity ARRAY PROCESSING | |
| CN103926572B (en) | A kind of clutter suppression method of battle array radar self adaptation subspace, airborne anon-normal side | |
| CN108375763B (en) | Frequency division positioning method applied to multi-sound-source environment | |
| CN108549059B (en) | Low-altitude target elevation angle estimation method under complex terrain condition | |
| CN104777453B (en) | The Beam Domain Time-Frequency Analysis Method of naval vessel line spectrum noise source positioning | |
| CN110196414B (en) | Antenna beam pointing method based on compensation of antenna directional diagram errors | |
| CN105335336B (en) | A kind of robust adaptive beamforming method of sensor array | |
| CN108828522A (en) | A kind of method of underwater vessel radiated noise measurement method using vertical array LCMV Wave beam forming | |
| CN114236464B (en) | Coherent interference resisting method for polarization sensitive array self-adaptive polarization adjustment | |
| CN108225536B (en) | Based on hydrophone amplitude and the self-alignment robust adaptive beamforming method of phase | |
| CN113281727B (en) | Output enhanced beam forming method and system based on horizontal line array | |
| CN105022050A (en) | Underwater-acoustic-channel discrete noise source suppression method of multi-sensor array | |
| CN105005038A (en) | Improved acoustic vector array coherent source DOA estimation algorithm | |
| Gebbie et al. | Passive localization of noise-producing targets using a compact volumetric array | |
| Gerstoft et al. | Adaptive beamforming of a towed array during a turn | |
| CN110196427A (en) | A kind of target location algorithm based on apart from orientation dictionary | |
| CN108957389A (en) | A kind of real number field multi channel signals method for estimating target azimuth | |
| CN110081964A (en) | Underwater sound source position and power spectrum combined estimation method based on sparse spectrum fitting | |
| CN109541526A (en) | A kind of ring array direction estimation method using matrixing | |
| CN111983599A (en) | Target two-dimensional DOA estimation method based on azimuth-pitch dictionary | |
| CN115826004B (en) | Three-star cooperative direct positioning method based on two-dimensional angle and time difference combination | |
| CN109814065B (en) | Beam forming method based on phase factor weighting | |
| CN113325401B (en) | Distortion towing linear array signal reconstruction method based on line spectrum phase difference deblurring | |
| CN115902849A (en) | Deep sea sound source depth estimation method based on beam output intensity resampling | |
| Swanson et al. | Small-aperture array processing for passive multi-target angle of arrival estimation |
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 |