CN104780127A - Multi-path channel estimation method based on time delay-Doppler R-L (Richardson-Lucy) deconvolution - Google Patents

Abstract

The invention discloses a multi-path channel estimation method based on time delay-Doppler R-L (Richardson-Lucy) deconvolution. The method comprises steps as follows: firstly, constructing a time delay-Doppler model about a waveguide multi-path environment, and deriving the expectation of a sampling mutual ambiguity function, that is, two-dimensional convolution of an emission signal self-ambiguity function and a channel spread function; then starting from the information theory, deriving an R-L iterative deconvolution algorithm according to the minimum Csiszar identification criterion; finally, performing two-dimensional deconvolution on a mean value of the sampling mutual ambiguity function to eliminate ambiguity introduced by an emission signal to obtain channel estimation with high-resolution time delay and Doppler parameters. Simulation examples indicate that the channel spread function can be reconstructed from the ambiguous sampling mutual ambiguity function with the R-L deconvolution algorithm; the resolution and precision of the R-L deconvolution algorithm are obviously higher than those of a traditional time delay-Doppler two-dimensional matched filtering method.

Description

A kind of delay-Doppler R-L uncoiling multi-path channel method of estimation

Technical field

The invention belongs to underwater acoustic channel to estimate, be specifically related to a kind of delay-Doppler R-L uncoiling multi-path channel method of estimation.

Background technology

Under water in the application such as sound communication, Ocean Acoustic Tomography and underwateracoustic detection, all estimate interested to underwater acoustic channel, because the handling property accurately directly having influence on whole system of channel estimating, this just requires the channel estimating of highly reliable resolution.

Oceanic waveguide be one has sea, under have seabed, the multi-path channel that there is the sound velocity gradient with change in depth centre, transmit by the time delay outside amount and doppler spread after this channel, Received signal strength is expressed as the duplicate weighted sum of the different delay-Doppler frequency shift transmitted.The two extended attribute of this delay-Doppler of ocean channel can represent with spread function, it and time varying channel impulse response function Fourier transform pairs each other.Thus channel estimating is by the channel impulse response functions from change when transmitting and Received signal strength estimates originally, is reduced to the estimation to channel expansion function.

Estimate parameters such as the time delay of channel expansion function and Doppler, traditional straightforward procedure is delay-Doppler matched filtering method, namely calculates and transmits and the mutual ambiguity function of sampling of Received signal strength.The resolving power of delay-Doppler matched filtering depends primarily on the main lobe width of signal from ambiguity function, has high time delay and Doppler's resolving power if transmitted, and as drawing pin Fuzzy degree function, the time delay that can obtain and Doppler estimate.

Mathematically, the Received signal strength of multi-path channel can be modeled as with the expectation of the mutual ambiguity function transmitted the bidimensional convolution transmitted from ambiguity function and channel expansion function, and thus the resolving power of delay-Doppler matched filtering method is by the restriction from ambiguity function characteristic.

Method for paying out has direct liftering method, but the method exists ill-conditioning problem.For solving ill-conditioning problem, people, under criterion of least squares instructs, successively propose multiple method for paying out, as the convex optimization method etc. of Wiener filtering method, regularization least square method, input and output nonnegativity restrictions.Criterion of least squares makes the square-error between true value and estimated value minimum, and this method is just best in Gaussian noise background, and there is noise scale-up problem, as Wiener filtering, and the problems such as regularization difficulty.

Summary of the invention

For tradition simply based on the multi-path channel method of estimation of delay-Doppler matched filtering, its resolving power is subject to the restriction of the main lobe width from ambiguity function, the present invention utilizes the nonnegativity restriction of channel expansion function, the minimum Csiszar in information theory is adopted to differentiate criterion, derive Richardson-Lucy (R-L) uncoiling algorithm, propose a kind of delay-Doppler R-L uncoiling multi-path channel method of estimation, its resolving power and precision significantly improve.

The present invention utilizes delay-Doppler R-L uncoiling, is achieved through the following technical solutions:

1) transmit with transmitting transducer;

2) hydrophone Received signal strength is used;

3) calculate transmit from ambiguity function;

4) Received signal strength and the mutual ambiguity function of sampling transmitted is calculated;

5) R-L algorithm is utilized to carry out bidimensional deconvolution to the mutual ambiguity function of sampling, until reach given condition to stop iteration;

6) obtain high-resolution multi-path channel to estimate, the horizontal stroke that each peak value is corresponding, ordinate are that the time delay τ in each path of channel and Doppler frequency shift ν estimate respectively.

Utilize a kind of delay-Doppler R-L uncoiling multi-path channel method of estimation of the present invention, comprise the steps:

(1) described step 2) signal that receives can carry out modeling, concrete modeling process is: point source transmits, Received signal strength transmits and the convolution of channel time varying impulse response, and be superimposed with noise, be formulated as r (t)=s (t) * h (τ, t)+n (t), wherein s (t), r (t), h (τ, t), n (t) represent respectively transmit, Received signal strength, channel time varying impulse response and noise, * represent one dimension convolution, parameter t, τ represent time and time delay respectively.

Channel expansion function is the Fourier transform of channel time varying impulse response, represents, can be formulated as H (τ, v)=∫ h (τ, t) e with H (τ, ν) ^{-j2 π vt}dt, wherein, j is imaginary unit's parameter, and ν represents Doppler frequency shift.Thus Received signal strength can be expressed as

r(t)＝∫∫H(τ,v)s(t-τ)e ^j2πvtdτdv+n(t) (1)

R (t), s (t) in formula, n (t) represents Received signal strength respectively, transmits and noise.

(2) if the bandwidth of signal is much smaller than centre frequency, such as ratio between two is less than 0.1; And the ratio of speed of related movement and the velocity of sound is much smaller than the inverse of signal period and bandwidth product, Doppler effect can be represented by Doppler frequency shift, and thus the Received signal strength of multi-path channel is expressed as

$r\left(t\right)=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ρ}}_{k}s(t-{\mathrm{\τ}}_k)e^{j2\mathrm{\π}{v}_{k}t}+n\left(t\right)---\left(2\right)$

In formula, k represents kth paths, and K represents total number of path; ρ _krepresent the amplitude-phase change of kth paths; τ _krepresent the time delay of kth paths; ν _krepresent the Doppler frequency shift of kth paths.

Formula (2) and formula (1) contrast learns that the channel expansion function under arrowband low-speed conditions is

$H(\mathrm{\τ},v)=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ρ}}_k\mathrm{\δ}(\mathrm{\τ}-{\mathrm{\τ}}_k,v-v_k)---\left(3\right)$

(3) described step 3) transmit be defined as follows from ambiguity function: the arrowband ambiguity function of finite energy signal s (t) is defined as

$\mathrm{\χ}(\mathrm{\τ},v)={\∫}_{-\∞}^{\∞}s\left(t\right)s^{*}(t-\mathrm{\τ})e^{-j2\mathrm{\πvt}}\mathrm{dt}---\left(4\right)$

In formula, subscript ^*represent conjugation, χ (τ, ν) is arrowband ambiguity function, weighs the similitude with the duplicate of its time delay, frequency displacement that transmits, and quantized signal is in the resolution capability of time delay and Doppler's bidimensional.

(4) described step 4) in Received signal strength and the mutual ambiguity function of sampling that transmits be defined as follows:

${\mathrm{\χ}}_c(\mathrm{\τ},v)={\∫}_{-\∞}^{\∞}r\left(t\right)s^{*}(t-\mathrm{\τ})e^{-j2\mathrm{\πvt}}\mathrm{dt}---\left(5\right)$

In formula, χ _c(τ, ν) is the mutual ambiguity function of sampling, weighs the similitude between transmitting after Received signal strength and time delay, frequency displacement.Received signal strength definition (1) is substituted into the mutual ambiguity function of sampling, draw the bidimensional convolution from ambiguity function that the expectation of mutual ambiguity function of sampling is channel expansion function and transmits

E[χ _c(τ,ν)]＝H(τ,ν)**χ(τ,ν) (6)

In formula, E [χ _c(τ, ν)] be the expectation of mutual ambiguity function of sampling, * * is bidimensional convolution.Always be positioned at this character of initial point from ambiguity function peak value learn according to what transmit, the corresponding time delay in position of mutual ambiguity function peak value and Doppler frequency shift, time delay and Doppler's resolving power are determined by the main lobe width from ambiguity function.

(5) described step 5) the criterion of R-L algorithm be: in the real function space of non-negative, wish to make Discrimination Functions minimum.Discriminating wherein refers to that Csiszar differentiates, is defined as

$L(p\left(x\right),q\left(x\right))={\∫}_{-\∞}^{\∞}p\left(x\right)\log \frac{p\left(x\right)}{q\left(x\right)}\mathrm{dx}-{\∫}_{-\∞}^{\∞}[p\left(x\right)-q\left(x\right)]\mathrm{dx}---\left(7\right)$

In formula, p (x), q (x) represent the real function of two non-negative, and L (p (x), q (x)) represents that the Csiszar between these two functions differentiates.With r (x), s (y), h (x|y), n (x) represents the output of filter, input, impulse response function and noise respectively, for one dimension convolution model r (x)=s (y) * h (x|y)+n (x), uncoiling task is in non-negative space, and known output and impulse response function ask input to estimate, the Csiszar met between output and output estimation differentiates minimum.If use represent that input is estimated, use represent output estimation, its computing formula is meet Csiszar is differentiated minimum input is estimated all must meet Kuhn-Tucker condition,

${\&Integral;}_{-\&infin;}^{\&infin;}\frac{h\left(x\right|y)}{\hat{r}\left(x\right)}r\left(x\right)\mathrm{dx}\left\{\begin{array}{c}=H_0\left(y\right),s\left(y\right)>0\\ <H_0\left(y\right),s\left(y\right)=0\end{array}\right.---\left(8\right)$

Wherein the analytic solutions meeting Kuhn-Tucker condition do not exist or are difficult to ask, and thus seek the iterative solution meeting this condition,

${\hat{s}}^{(i+1)}\left(y\right)=\frac{{\hat{s}}^{\left(i\right)}\left(y\right)}{H_0\left(y\right)}{\&Integral;}_{-\&infin;}^{\&infin;}\frac{h\left(x\right|y)}{{\&Integral;}_{-\&infin;}^{\&infin;}h\left(x\right|y){\hat{s}}^{\left(i\right)}\left(y\right)\mathrm{dy}}r\left(x\right)\mathrm{dx}---\left(9\right)$

In formula, i represents iterations, with represent the input estimation of iteration i+1 time and i time respectively.Proved that R-L algorithm is iteration convergence, and R-L algorithm is equally applicable to bidimensional uncoiling.Known transmit from ambiguity function, according to Received signal strength and transmit calculate mutual ambiguity function of sampling average after, obtain that the estimation of channel expansion function is iterative is

${\hat{H}}^{(i+1)}(\mathrm{\&tau;},v)={\hat{H}}^{\left(i\right)}(\mathrm{\&tau;},v)[\frac{\mathrm{\&chi;}(\mathrm{\&tau;},v)}{\mathrm{\&chi;}(\mathrm{\&tau;},v)**{\hat{H}}^{\left(i\right)}(\mathrm{\&tau;},v)}**E[{\mathrm{\&chi;}}_c(-\mathrm{\&tau;},-v)\left]\right]---\left(10\right)$

In formula, with represent the estimation of iteration i+1 time and the channel expansion function of i time respectively, this bidimensional solution rolling eliminates normalization coefficient ∫ ∫ χ (τ, ν) d τ d ν.When reaching stopping criterion for iteration, stop iteration uncoiling, the horizontal stroke that in the estimated result obtained, each peak value is corresponding, vertical two-dimensional position are that the time delay τ in each path of channel and Doppler frequency shift v estimate.

The present invention by the mutual ambiguity function of sampling is carried out bidimensional uncoiling eliminate transmit fuzzyly estimate channel expansion function, the delay-Doppler uncoiling algorithm proposed does not require to transmit and has desirable ambiguity function, but the time delay progressively realized by the method for iteration and Doppler are estimated.

Accompanying drawing explanation

Fig. 1 represents the workflow diagram that a kind of delay-Doppler R-L of the present invention uncoiling multi-path channel is estimated;

Fig. 2 represent one transmit from ambiguity function;

Fig. 3 represents that one transmits and the mutual ambiguity function of Received signal strength, wherein

Fig. 3 (a) represents the 3-D view of mutual ambiguity function,

Fig. 3 (b) represents time delay and the Doppler frequency shift two-dimensional view of mutual ambiguity function;

Fig. 4 represents the estimation of the channel expansion function utilizing R-L uncoiling to obtain, wherein

Fig. 4 (a) represents the 3-D view of channel expansion Function Estimation,

Fig. 4 (b) represents time delay and the Doppler frequency shift two-dimensional view of channel expansion Function Estimation;

Fig. 5 represents the relation between uncoiling mean square error and iterations.

Embodiment

Below in conjunction with accompanying drawing and instantiation, the present invention will be further described, to verify validity of the present invention.Fig. 1 represents the workflow diagram that a kind of delay-Doppler R-L of the present invention uncoiling multi-path channel is estimated, concrete implementation process is as follows:

(1) select to transmit, as adopted QPSK signal, symbol rate is 1KHz, and for estimating that the symbol numbers of channel is 200, calculate it from ambiguity function, its two-dimensional view as shown in Figure 2.The resolving power of delay-Doppler matched filtering be inversely proportional to from the main lobe width of ambiguity function, specifically, time delay resolving power and signal bandwidth are inversely proportional to, and Doppler's resolving power and signal pulsewidth are inversely proportional to.

(2) set channel parameter, number of path is set to 4, and design parameter is as table 1.The velocity of sound is 1500m/s, and introduce white Gaussian noise, signal to noise ratio is 5dB, parameter is substituted into formula (2) and obtains Received signal strength.

Table 1 path parameter

(3) calculate the mutual ambiguity function transmitted with Received signal strength, as shown in Figure 3, (a) (b) represents 3-D view and two-dimensional view respectively, represents the position in each path in table 1 in figure (b) by black circles.As shown in Figure 3 signal be about 1ms from the time delay resolving power of ambiguity function, path 3 and 4 can be differentiated, path 1 and 2 can not be differentiated; Doppler's resolving power is about 5Hz, can not estimate the Doppler of 4 paths well.Simultaneously because the Doppler frequency shift interval in path 1 and 2 only differs 0.5Hz, and time delay spacing difference 1ms, the main lobe at two paths places has superposed mutually, can not tell path 1 and 2 completely.

(4) the sampling mutual ambiguity function of R-L algorithm iteration formula (10) to above-mentioned steps (3) is utilized to carry out bidimensional uncoiling, the channel expansion Function Estimation that iteration is 100 times as shown in Figure 4, a () (b) represents 3-D view and two-dimensional view respectively, Fig. 4 (b) black circles represents the position in each path in table 1.Uncoiling result can tell four paths completely as shown in Figure 4, and it is consistent that location parameter and table 1 set, and illustrates that R-L uncoiling algorithm can be estimated channel expansion function well.Contrast with the delay-Doppler matched filtering method of mutual ambiguity function of sampling that directly adopts shown in Fig. 3, be all significantly improved based on the resolving power of the delay-Doppler method for paying out of R-L algorithm and precision.

(5) for further illustrating the Iterations of Multi of R-L algorithm, the graph of a relation between mean square error (Mean Square Error, MSE) and iterations is drawn, as Fig. 5.MSE is defined as the mean square error between estimates of parameters and parameter true value, and estimates of parameters is uncoiling result, and be amplitude attenuation value according to table 1 in respective path place parameter true value, other positions are zero.Can find the increase along with iterations in figure, MSE dullness reduces.

Claims (6)

1. a delay-Doppler R-L uncoiling multi-path channel method of estimation, is characterized in that, comprising:

1) transmit with transmitting transducer;

2) hydrophone Received signal strength is used;

3) calculate transmit from ambiguity function;

4) Received signal strength and the mutual ambiguity function of sampling transmitted is calculated;

5) Richardson-Lucy (R-L) algorithm is utilized to carry out bidimensional deconvolution, until the iteration ends that satisfies condition to the mutual ambiguity function of sampling;

6) obtain multi-path channel by iteration to estimate, the horizontal stroke that each peak value is corresponding, ordinate are that the time delay τ in each path of channel and Doppler frequency shift ν estimate respectively.

2. multi-path channel method of estimation as claimed in claim 1, it is characterized in that, described Received signal strength is:

r(t)＝∫∫H(τ,v)s(t-τ)e ^j2πvtdτdv+n(t)

H(τ,v)＝∫h(τ,t)e ^-j2πvtdt

In formula, j is imaginary unit, and t represents the time, τ is time delay, and r (t) is Received signal strength, and s (t) is for transmitting, s (t-τ) is the signal after time delay τ, and n (t) is noise, and ν is Doppler frequency shift, H (τ, ν) is channel expansion function, h (τ, t) be channel time varying impulse response, h (τ, t) and H (τ, ν) is a pair fourier transform pair.

3. multi-path channel method of estimation as claimed in claim 1, is characterized in that, described in transmit from ambiguity function, adopt Woodward arrowband ambiguity function, be defined as:

$\mathrm{\&chi;}(\mathrm{\&tau;},v)={\&Integral;}_{-\&infin;}^{\&infin;}s\left(t\right)s^{*}(t-\mathrm{\&tau;})e^{-j2\mathrm{\&pi;vt}}\mathrm{dt}$

In formula, χ (τ, ν) is arrowband ambiguity function, subscript ^*represent conjugation.

4. multi-path channel method of estimation as claimed in claim 1, it is characterized in that, the mutual ambiguity function of described sampling is:

${\mathrm{\&chi;}}_c(\mathrm{\&tau;},v)={\&Integral;}_{-\&infin;}^{\&infin;}r\left(t\right)s^{*}(t-\mathrm{\&tau;})e^{-j2\mathrm{\&pi;vt}}\mathrm{dt}$

In formula, χ _c(τ, ν) is the mutual ambiguity function of sampling.

5. multi-path channel method of estimation as claimed in claim 4, it is characterized in that, the mutual ambiguity function of sampling described in being substituted into by Received signal strength, draws the bidimensional convolution from ambiguity function that the expectation of mutual ambiguity function of sampling is channel expansion function and transmits

E[χ _c(τ,ν)]＝H(τ,ν)**χ(τ,ν)

In formula, E [χ _c(τ, ν)] be the expectation of mutual ambiguity function of sampling, H (τ, ν) is channel expansion function, and * * represents bidimensional convolution, χ (τ, ν) be transmit from ambiguity function.

6. multi-path channel method of estimation as claimed in claim 5, it is characterized in that, the estimation of the channel expansion function obtained by delay-Doppler bidimensional R-L iteration uncoiling algorithm meets relational expression:

${\hat{H}}^{(i+1)}(\mathrm{\&tau;},v)={\hat{H}}^{\left(i\right)}(\mathrm{\&tau;},v)[\frac{\mathrm{\&chi;}(\mathrm{\&tau;},v)}{\mathrm{\&chi;}(\mathrm{\&tau;},v)**{\hat{H}}^{\left(i\right)}(\mathrm{\&tau;},v)}**E[{\mathrm{\&chi;}}_c(-\mathrm{\&tau;},-v)\left]\right]$

In formula, i represents iterations, with represent the estimation of iteration i+1 time and the channel expansion function obtained for i time respectively.