CN106603036A - Adaptive time delay estimation method based on low-order interpolation filter - Google Patents

Abstract

The invention provides an adaptive time delay estimation method based on a low-order interpolation filter, which can synthesize any continuously-variable non-integral period sampling delay in real-time and optimal mode. The any k time is set as the initial estimation time, and the initial value of the time delay estimation value is set according to prior information. The time delay estimation value is decomposed to obtain an integral-period delay component and a decimal-period delay component. An index position variable is calculated. The interpolation filter is adopted to carry out interpolation operation on observation signals, and estimation signals are generated. According to the observation signals and the generated estimation signals, the real-time error, the lag error and the lead error at the any kth time are calculated. The instantaneous gradient of an LMS adaptive algorithm objective function is calculated, and finite difference is used to approximate partial derivative operation. The time delay estimation value at the any (k+1)th time is calculated, the above steps are iterated, judgment is carried out after each time of iteration, and when the objective function is smaller than a set constant, the time delay estimation value at the iteration at the time serves as the final estimation result for the time-varying delay.

Description

A kind of novel adaptive time delay estimation based on low order interpolation filter

Technical field

The invention belongs to signal time delay estimation technique field, and in particular to a kind of self adaptation based on low order interpolation filter Delay time estimation method.

Background technology

Time delay is estimated as the theoretical important component part of modern signal processing, is always compared in digital processing field Relatively active research direction, is widely used in the technologies such as radar, sonar, navigation, remote measurement, global location, seismology, biomedicine Field.Delay time estimation method species is various, mainly including method of correlation, lowest mean square (LMS-Least Mean Square) method, pole Maximum-likelihood method, feature structure method, High order statistics etc..In above method, the novel adaptive time delay estimation based on LMS due to Without the need for the prior information of known noise and signal, have the advantages that robustness is good, amount of calculation is little, can dynamic tracking time delay change, The extensive concern for just receiving Chinese scholars once coming out.

In LMS novel adaptive time delay estimations, ETDE (Explicit Time Delay Estimator) method is recognized For be most classical, one of the method that is most widely used.The core concept of ETDE is that time delay is modeled as into the little of Sinc type Number time delay FIR filter, so as to directly be updated to obtain the non-integral multiple sampling period to time delay in adaptive algorithm Delay Estima-tion value.Compared to traditional LMSTDE (Least Mean Square Time Delay Estimator), ETDE To the interpolating operations of wave filter weight coefficient in negligible LMSTDE, dash forward with estimated accuracy is high, amount of calculation is little and real-time is good etc. Go out advantage.However, ETDE has proved to be a kind of biased estimation under finite filter length, time delay estimated bias can be with letter Make an uproar and increase than the reduction with filter order.By increasing an adaptive gain control, ETDGE in ETDE (Explicit Time Delay and Gain Estimator) can realize the unbiased of time delay under broader filter length Estimate.Despite this, ETDE and ETDGE algorithms can only also be adapted well to the time delay between Whole frequency band noise class signal estimating. For the time delay of band-limited signal is estimated, because Sinc type decimal filtering wave by prolonging time utensils have sizable passband ripple, therefore The Time delay Estimation Accuracy of ETDE and ETDGE algorithms is by severe exacerbation.During in view of the LMS self adaptations for being based on decimal filtering wave by prolonging time device The estimated accuracy for prolonging method of estimation depends primarily on the pass-band performance of decimal filtering wave by prolonging time device, is that this Chinese scholars is generally adopted Lagrange interpolations replace Sinc interpolations and estimate carrying out the time delay between band-limited signal, such as based on the ETDE of Lagrange interpolations, ETDE based on AMDF (Average Magnitude Difference Function) etc..Comparatively speaking, Lagrange types are little Number filtering wave by prolonging time device has maximally-flat frequency characteristic at low frequency, approaches when filter order is identical and need to only carry out low pass When, its delay precision is much better than Sinc type decimal filtering wave by prolonging time devices.If the mid frequency of band-limited signal is, it is known that generally can be by decimal The modulation of filtering wave by prolonging time device is arrived at signal center frequency, with relatively low filter order, around the limited of signal center frequency Higher delay precision is provided in bandwidth.At present, mainly included based on the LMS delay time estimation methods of modulation decimal filtering wave by prolonging time device METDE (Modulated ETDE), MLETDE (Modulated Lagrange ETDE) and MMLETDE (Mixed Modulated Lagrange ETDE) etc..

However, said method is mostly based on single decimal filtering wave by prolonging time device (or interpolater), can only preferably be applied to The time-varying delays of a small range change are estimated.For large-scale time-varying delays are estimated, it is inevitably present as follows Two aspect problems：On the one hand, the time delay of existing method is estimated and following range directly depends on the exponent number of decimal filtering wave by prolonging time device, Estimate that it must be introduced into the decimal filtering wave by prolonging time device of high-order, and this undoubtedly considerably increases fortune to adapt to large-scale time-varying delays The complexity of calculation；On the other hand, the optimization time delay interval of decimal filtering wave by prolonging time device is only located at the center of filter impulse response Place, although can increase time delay estimation and following range by increasing filter order, but when time delay change exceedes this optimization area Between when, Time delay Estimation Accuracy is by severe exacerbation.

The content of the invention

In view of this, the invention provides a kind of novel adaptive time delay estimation based on low order interpolation filter, can In real time, the variable non-integer-period sampled time delay of arbitrary continuation is most optimally synthesized.

In order to achieve the above object, the technical scheme is that：During a kind of self adaptation based on low order interpolation filter Prolong method of estimation, source signal is s (k), and the method is used in measurement noise n₁(k) and n₂In the presence of (k), connect according to sensor The change of time-varying delays D (k) is estimated and tracked to two observation signals x (k) and y (k) for receiving；Wherein k is time variable.

Step 1：If using any k moment as the estimated initial moment, the estimated value of time-varying delays D (k) isAccording to elder generation Test information settingInitial value be

Step 2：It is rightDecomposed：Obtain delay component complete cycleWith Fractional delay componentI.e.Wherein,Downward floor operation is represented,For integer,For decimal and presence

Step 3：According toComputation index location variableI.e.

Step 4：According toWith the mid frequency ω of source signal s (k)₀, using interpolation filter to observation signal x K () carries out interpolative operation, produceThe interpolation formula of the interpolation filter is specific as follows：

N is the exponent number of the interpolation filter, that is, be by sample index location variableIt is determined that any k moment interpolation when Quantity N of required observation signal sampled point, c (m, p) is the complex value coefficient of setting；M be interpolation variable, m ∈ [- M₁,M₂],p∈ [0, N-1] and there is N=M₁+M₂；

Step 5：Calculate what is generated according to observation signal y (k) and step 4Calculate any kth moment Immediate errorI.e.：

Step 6：Taking lagged value isIt is with advance valueH be default step-length because Son；It is rightWithBy the way of step 2～step 5, hysteresis error is respectively obtainedWith advanced error

According to immediate errorHysteresis errorWith advanced errorCalculate LMS adaptive Answer the temporary gradients of algorithm object functionWherein derivative operation is approached using finite difference.

Step 7：According to the time delay estimated value at kth momentWith the temporary gradients of object functionPress Formula iterates to calculate the time delay estimated value at any moment of kth+1I.e.μ is convergence The factor.

Step 8：WillAs in step 2Repeat step 2 is iterated to step 7, sentences after each iteration It is disconnected, object functionWhether less than setting constant ε, if then in this iterationIt is designated as time-varying delays D (k) Final estimated result, the time delay for otherwise calculating subsequent time by step 7 is estimated to and as in step 2Continuation changes Generation.

Further, c (m, p) is the complex value coefficient of setting in step 4, and establishing method is as follows：

Wherein adoptObtain c₀The value of (m, p).

Further, in step 6, the temporary gradients of LMS adaptive algorithm object functionsSpecially：

Wherein,For object function, Re { } is represented and is taken real operation, and * represents complex conjugate operation.

Beneficial effect：

Interpolation filter in the present invention has following features：1st, for arbitrarily large non-integer-period time delay, interpolation Exponent number is unrelated with time delay size；2nd, interpolation coefficient is only dependent upon the decimal delay component of any non-integer-period time delay, and decimal prolongs When part be always positioned at digit pulse response center；3rd, the complex value coefficient adopted in interpolation formula of the invention is fixed not Become, and without the need for calculating in real time.Understood based on These characteristics, the present invention puies forward the method for online interpolation and can adopt relatively low interpolation rank Several interpolation filters in real time, most optimally synthesizes the variable non-integer-period sampled time delay of arbitrary continuation realizing, and compares In traditional ETDE methods based on decimal filtering wave by prolonging time device, under identical filter order, institute's extracting method of the present invention can be carried For bigger Delay Estima-tion and following range, and higher Time delay Estimation Accuracy.

Description of the drawings

Fig. 1 is the theory diagram of the present invention.

Specific embodiment

Develop simultaneously below in conjunction with the accompanying drawings embodiment, describes the present invention.

Embodiment 1, the inventive method is achieved through the following technical solutions：

It is a kind of based on low order interpolation filter on a large scale, high-accuracy self-adaptation time-varying delays method of estimation, it is substantially real Apply process as follows：

If the mathematical model that two sensors receive discrete observation signal is represented by

Wherein, k is discrete-time variable, and s (k) is limited bandwidth, mid frequency is ω₀Source signal, n₁(k) and n₂(k) For zero-mean, orthogonal additive white Gaussian, x (k) and y (k) is observation signal, and D (k) is time-varying delays on a large scale.

Institute's extracting method of the present invention is sought in measurement noise n₁(k) and n₂In the presence of (k), received according to sensor Two observation signals x (k) and y (k) in real time, accurately estimate and track the change of time delay D (k).It is specific as follows：

Step 1：If using any kth moment as the estimated initial moment, being set according to prior informationAs (big model Enclose) the estimation initial value of time-varying delays D (k)

Step 2：It is rightDecomposed as the following formula, be calculated delay component complete cycleWith Fractional delay componentIt is specific as follows：

I.e.

Wherein,Downward floor operation is represented,For integer,For decimal and presence

Step 3：According to delay component complete cycleComputation index location variableI.e.

Step 4：According to index position variableDecimal time delayWith the mid frequency ω of source signal s (k)₀, by (5) formula Interpolative operation is carried out to observation signal x (k), is producedIt is specific as follows：

Wherein, c (m, p) is fixed multiple interpolation filter coefficient, and computational methods are as follows：

Wherein, m ∈ [- M₁,M₂], p ∈ [0, N-1], and there is N=M₁+M₂

The step principle is as follows：

The purpose of this step be observation signal x (k) is carried out on a large scale, the non-integer-period time delay of high accuracy, continuous variableIt is calculatedIt is arbitrarily band-limited continuous according to knowable to the sample reconstruction of band-limited signal is theoretical Time signal x (t) can be by its sample sequence { x (nT_s) be reconstructed as the following formula (without loss of generality, it is assumed that T_s=1)：

Wherein, h_IT () is the shock response of ideal interpolation wave filter continuous time, with following form：

Continuous time signal x (t) is carried outTime delay, then with t=kT_s=k carries out resampling, can obtain

It can be seen that, time delayed signalThe time delay that can be responded by the sample sequence of observation signal and ideal interpolation filter impulse Sampling is obtained by convolution algorithm.However, formula (11) is related to unlimited item summation operation, engineering cannot be realized.For this purpose, pushing away below The one kind for leading (11) formula approaches form.

It is rightDecomposed as the following formula

Wherein,For delay component complete cycle,For Fractional delay component.

Index of definition location variable is

Interpolation filter variable is

(12)-(14) formula is substituted into into (11) formula, can be obtained

Wherein,To approach interpolation filter, index position variableFor determining any kth time interpolation when institute The N number of observation signal sampled point for needing, decimal time delayDigit pulse for calculating interpolation filter is responded

The present invention determines the digit pulse response of interpolation filter with reference to Lagrange interpolationsSimultaneously in order to just In the continuous control of time delay, using decimal time delayN rank multinomials pairApproached, i.e.,

Wherein, c₀(m, p) is fixed real-valued Coefficients of Approximation, with decimal time delayIt is unrelated.

In order to ensure interpolation filter in source signal mid frequency ω₀Higher time delay essence is provided in neighbouring limited-bandwidth Degree, can be modulated as the following formula, i.e., to (16) formula

Wherein,For the interpolation filter coefficient of modulation, c (m, p) is fixed complex value Coefficients of Approximation, and is existed

Wherein, m ∈ [- M₁,M₂],p∈[0,N-1]。

By the interpolation filter coefficient of modulationIn replacing (15) formulaFinal interior of this step can be obtained Slotting formula is as follows：

It can be seen from above formula, for arbitrarily large non-integer-period time delayFor, the present invention carries interpolation formula to be had Following distinguishing feature：

1st, the interpolation filter of modulationWithUnrelated, interpolation filter exponent number N can not be withIncrease And increase；

2nd, the digit pulse response of interpolation filter is only dependent uponDecimal delay componentIt is non-whole for arbitrary Cycle time delayDecimal time delayDecimal filtering wave by prolonging time device is located at foreverOptimization time delay it is interval；

3rd, complex value coefficient c (m, p) immobilizes, and without the need for calculating in real time, facilitates implementation the continuous control of time delay.

Understood based on These characteristics, the present invention carries interpolation formula can be under relatively low interpolation filter exponent number, in real time, most Optimally synthesize the variable non-integer-period sampled time delay of arbitrary continuation.

Step 5：Calculate what is generated according to observation signal y (k) and step 4Calculate any kth moment Immediate errorI.e.

Step 6：While step 2 to step 5 is carried out, make respectivelyWithReplace In step 2By step 2 to step 5 method, while being calculated hysteresis errorWith advanced errorI.e.

Wherein, h is step factor.

According to immediate errorHysteresis errorWith advanced errorCalculate LMS adaptive Answer the temporary gradients of algorithm object functionI.e.

Wherein,For object function, Re { } is represented and is taken real operation, and " * " represents complex conjugate operation.

The step principle is as follows：

In LMS adaptive algorithms, instantaneous mean square error is generally taken as object function, i.e.,

Wherein,For object function, " * " represents complex conjugate operation.

So, the temporary gradients of object functionCan be calculated as follows：

Wherein, symbolTable takes derivative operation.

In order to reduce the operand of adaptive iteration algorithm, the present invention approaches derivative operation using finite difference, i.e.,

And exist

(25) formula is substituted into into (24) formula to understand, the temporary gradients of object functionCan approximate representation be

Step 7：According to the time delay estimated value at kth momentWith the temporary gradients of object functionPress Formula iterates to calculate the time delay estimated value at any moment of kth+1I.e.

Wherein, μ is convergence factor, for controlling the stability and convergence rate of adaptive algorithm.

Step 8：WillAs in step 2Repeat step 2 is iterated to step 7, sentences after each iteration It is disconnected, object functionWhether less than setting constant ε, if then in this iterationIt is designated as time-varying delays D (k) Final estimated result, the time delay for otherwise calculating subsequent time by step 7 is estimated to and as in step 2Continuation changes Generation.

It can be seen that, through step 1 to step 8, the present invention is in measurement noise n₁(k) and n₂In the presence of (k), according to observation letter Number x (k) and y (k) in real time, accurately estimate and have tracked the change of time delay D (k).

To sum up, presently preferred embodiments of the present invention is these are only, is not intended to limit protection scope of the present invention.It is all Within the spirit and principles in the present invention, any modification, equivalent substitution and improvements made etc. should be included in the protection of the present invention Within the scope of.

Claims (3)

1. a kind of novel adaptive time delay estimation based on low order interpolation filter, it is characterised in that source signal is s (k), the party Method is used in measurement noise n₁(k) and n₂In the presence of (k), two observation signals x (k) received according to sensor and y (k) Estimate and track the change of time-varying delays D (k)；Wherein k is time variable；

Step 1：If using any k moment as the estimated initial moment, the estimated value of time-varying delays D (k) isBelieved according to priori Breath settingInitial value be

Step 2：It is rightDecomposed：Obtain delay component complete cycleWith Fractional delay componentI.e.Wherein,Downward floor operation is represented,For integer,For decimal and presence

Step 3：According toComputation index location variableI.e.

Step 4：According toWith the mid frequency ω of source signal s (k)₀, observation signal x (k) is entered using interpolation filter Row interpolative operation, producesThe interpolation formula of the interpolation filter is specific as follows：

$\begin{array}{c}\hat{y}\left(k\right)=x\left(k-\hat{D}\left(k\right)\right)\\ \≅\underset{p=0}{\overset{N-1}{\mathrm{\Σ}}}\[\underset{m=-M_1}{\overset{M_2}{\mathrm{\Σ}}}c\left(m,p\right)x\left({\hat{m}}_k-m\right)\]{\hat{d}}_k^p{e}^{-{\mathrm{j\ω}}_0{\hat{d}}_k}\end{array};$

N is the exponent number of the interpolation filter, that is, be by sample index location variableIt is determined that any k moment interpolation when it is required Observation signal sampled point quantity N, c (m, p) be setting complex value coefficient；M be interpolation variable, m ∈ [- M₁,M₂],p∈[0, N-1] and there is N=M₁+M₂；

Step 5：Calculate what is generated according to observation signal y (k) and step 4Calculate any kth moment i.e. Time errorI.e.：

Step 6：Taking lagged value isIt is with advance valueH is default step factor；It is rightWithBy the way of step 2～step 5, hysteresis error is respectively obtainedWith advanced error

According to immediate errorHysteresis errorWith advanced errorCalculate LMS self adaptations to calculate The temporary gradients of method object functionWherein derivative operation is approached using finite difference；

Step 7：According to the time delay estimated value at kth momentWith the temporary gradients of object functionIteration as the following formula Calculate the time delay estimated value at any moment of kth+1I.e.μ is convergence factor；

Step 8：WillAs in step 2Repeat step 2 is iterated to step 7, judges after each iteration, target FunctionWhether less than setting constant ε, if then in this iterationIt is designated as finally estimating for time-varying delays D (k) Meter result, the time delay for otherwise calculating subsequent time by step 7 is estimated to and as in step 2Continue iteration.

2. a kind of novel adaptive time delay estimation based on low order interpolation filter as claimed in claim 1, it is characterised in that C (m, p) is the complex value coefficient of setting in the step 4, and establishing method is as follows：

$c(m,p)=c_0(m,p)e^{{\mathrm{j\ω}}_{0}m};$

Wherein adoptObtain c₀The value of (m, p).

3. a kind of novel adaptive time delay estimation based on low order interpolation filter as claimed in claim 1, it is characterised in that In the step 6, the temporary gradients of LMS adaptive algorithm object functionsSpecially：

Wherein,For object function, Re { } is represented and is taken real operation, and * represents complex conjugate operation.