Subspace-based blind identification algorithm of ocean underwater acoustic channel for multi-channel FIR filter
11050492 · 2021-06-29
Assignee
Inventors
Cpc classification
H04B11/00
ELECTRICITY
International classification
Abstract
The disclosure provides a subspace-based blind identification algorithm of an ocean underwater acoustic channel for multi-channel fir filter, which adopts a technical solution that a channel impulse response coefficient is calculated by quadratic minimization. The disclosure has beneficial effects that estimation precision can be met when using a proper number of samples, and especially when a few noise vectors are used for estimating channel parameters, so that calculation amount is greatly reduced.
Claims
1. A subspace-based blind identification algorithm of an ocean underwater acoustic channel for multi-channel FIR filter, comprising: calculating a channel impulse response coefficient by quadratic minimization according to a formula (16):
{tilde over (q)}(h)=h.sup.H{tilde over (Q)}h (17); solving a maximum eigenvalue of {tilde over (Q)} which makes the formula (17) maximized; estimating a correlation matrix R.sub.x by replacing an overall average of the correlation matrix of communication signals with a sample estimation due to the overall average of the correlation matrix of communication signals being unknown in actual communication, wherein
R.sub.x=Sdiag(λ.sub.0,λ.sub.1, . . . λ.sub.L+K−1)S.sup.H+σ.sup.2GG.sup.H (8); wherein, columns of matrix S are expanded into signal subspaces, and columns of matrix G are expanded into noise subspaces; and a calculation of a quadratic form: solving min[q(h)]=min[h.sup.TQh] according to formula (15), and solving min(Q.sub.i)=min(ĝ,ĝ.sub.i.sup.H), wherein h=Q.sup.−1c.
Description
BRIEF DESCRIPTION OF THE DRAWINGS
(1)
(2)
(3)
(4)
(5)
(6)
(7)
DETAILED DESCRIPTION
(8) Hereinafter, the disclosure will be described in detail with reference to specific embodiments.
(9) 1. Blind Identification Model of Multi-Channel FIR Filter
(10) It is assumed that discrete signals transmitted by a signal source are s(k), signals observed after s(k) have been transmitted through channels, filtered and processed are defined as x(k), an impulse response of the channel is h(k) including comprehensive characteristics of transmission, reception, processing of channel and the like, h(k) is linear, and a noise n(k) superimposed during transmission is a band-limited stationary process.
(11) The subspace method relies on two assumptions that: (1) the channel impulse response h(k) has a finite support, i.e. the channel can be represented by an FIR filter; (2) multiple measurements are possible within a character transmission period T.
(12) There are several methods to implement multiple measurements, one of which is adopting multiple sensors to sample signals received by the sensors with the period T, another of which is oversampling a signal received by a single sensor, or a comprehensive application of both. Both the method adopting multiple sensors and the oversampling method adopting a single sensor are equivalent, and they are equivalent to a single-input multiple-output system, as shown in
(13) Input signal sequence {s(k)} is output through different linear channels {h.sub.m(k)}, m=1, 2, . . . , M. An output of any channel can be represented as:
(14)
(15) wherein “.Math.” represents convolution, H(k)=[h.sub.1(k), h.sub.2(k), . . . , h.sub.M (k)].sup.T, X(k)=[x.sub.1(k), x.sub.2(k), . . . , x.sub.M(k)].sup.T, and N(k)=[n.sub.1(k), n.sub.2(k), . . . , n.sub.M(k)].sup.T.
(16) If the SIMO system is a finite impulse response system (FIR) with a length of L, when an observation length is K, an output of the SIMO system is:
X.sub.K=H.sub.KS.sub.K+N.sub.K (2);
(17) wherein each matrix in formula (2) is expressed as:
(18)
(19) wherein X.sub.k and N.sub.k are both K×M observation signal matrices, H.sub.K is an M×K×(K+L) filter matrix, and S.sub.k is a (K+L)×1 vector.
(20) It has been proved in Document [71] (i.e. [71] Li Y, Ding Z. Blind Channel Identification Based on Second-order Cyclosttationary statistics. ICASSP, 1993, 4:81-84.) that channel parameters of the SIMO system can be identified by using the second-order statistics output by the system, when M channels meet two following conditions: (1) x.sub.m.sub.
(21) 2. Subspace-Based Blind Identification of Multi-Channel FIR Filter
(22) The blind identification is a process of estimating the channel impulse response coefficients H.sub.k only according to the observed data X.sub.k. The formula (2) can be rewritten as:
X.sub.k=H.sub.kS.sub.k+N.sub.k (3).
(23) The blind identification is based on an MK×MK autocorrelation matrix R.sub.x according to Tong's method (see Document [72] L. Tong, G. Xu, T. Kailath. Blind identification and equalization based on second-oreder statistics: a time domain approach. IEEE Trans. Inform. Theory, 1994, 40(3): 340-350.), i. e:
R.sub.x=E(X.sub.kX.sub.k.sup.H) (4).
(24) Formula (5) can be built when the formula (3) is substituted into the formula (4) due to an assumption of additive noise being independent from transmitted signal sequence:
R.sub.x=H.sub.kR.sub.sH.sub.k.sup.H+R.sub.N (5),
(25) wherein R.sub.s=E(S.sub.kS.sub.k) and R.sub.N=E(N.sub.kN.sub.k.sup.H) represents an autocorrelation matrix of a transmitted signal vector S.sub.k and that of a noise vector N.sub.k respectively. The autocorrelation matrix R.sub.5 of the signal source has dimensions of (K+L)×(K+L) and is assumed to be full-rank but unknown. The noise autocorrelation matrix has dimensions of KM×KM and is assumed as R.sub.N=∝.sup.2I for simplicity, wherein a represents a variance of N.sub.k, I represents a KM×KM unit matrix.
(26) λ.sub.0≥λ.sub.1≥ . . . λ.sub.KM−1 is assumed to be eigenvalues of R.sub.x. A signal part of R.sub.x, i.e. H.sub.kH.sub.sH.sub.k.sup.H, has a rank of K+L due to R.sub.s being full-rank, that is:
(27)
(28) Unit eigenvectors corresponding to the eigenvalues λ.sub.0≥λ.sub.1≥ . . . λ.sub.K+M−1 are denoted as S.sub.0, S.sub.1, . . . , S.sub.K+M−1, and unit eigenvectors corresponding to eigenvalues λ.sub.L+K, λ.sub.L+K+1, . . . λ.sub.MK−1 are denoted as G.sub.0, G.sub.1, . . . , G.sub.KM−(L+K)−1, and there is:
(29)
(30) wherein S has dimensions of MK×(L+K) and G has dimensions of MK×(MK−L−K).
(31) The matrix R.sub.s can then be represented as:
R.sub.x=Sdiag(λ.sub.0,λ.sub.1, . . . λ.sub.L+K−1)S.sup.H+σ.sup.2GG.sup.H (8);
(32) Wherein columns of matrix S are expanded into signal subspaces and columns of matrix G are expanded into noise subspaces. The noise subspace is an orthocomplement of the signal subspace, and the signal subspace is also a linear space expanded from columns of a filter H.sub.K. It can be known from the orthogonality of the noise and signal subspaces that the columns of the filter matrix H.sub.K are orthogonal to any vector within the noise subspaces, which can be represented as:
G.sub.i.sup.HH.sub.K=O,0≤i≤MK−L−K (9).
(33) Although the formula (9) is a linear equation with unknown coefficients, it provides a simple method for blind channel identification. The impulse response coefficients of channels can be uniquely determined by the noise subspace under a condition of at most one constant gain difference with appropriate conditions in the following theorem 1.
(34) Theorem 1: assuming that K≥L, a matrix H.sub.K−1 is column full rank (i.e., rank(H.sub.K−1)=L+K−1), and H′.sub.K is a non-zero filter matrix with same dimensions as H.sub.K, a necessary and sufficient condition that a numerical range of H′.sub.K is included within a numerical range of H.sub.K is that H′.sub.K is proportional to H.sub.K.
(35) The theorem 1 provides a method for estimating channel parameters, wherein the identification of channels depends on a particular structure of H.sub.K under a condition that R.sub.s is unknown and full-rank.
(36) 3. Subspace-Based Parameter Estimation
(37) Orthogonal condition formula (9) is a linear function of channel coefficient, but in practice only a sample estimation Ĝ.sub.i of the noise eigenvectors can be used for parameter estimation, so that formula (9) should be solved only in the least square sense. A quadratic cost function is defined as below:
(38)
(39) It is easily derived from the following lemma 1 that q(h) depends on h rather than H.sub.K. To this end, the following symbols are introduced, and the structures of h and H.sub.K are the same therewith.
(40) It is assumed that v.sup.(0), v.sup.(1), . . . , v.sup.(M−1) are M arbitrary K×1 vectors, denoted as MK×1 vectors of v=[v.sup.(0), v.sup.(1), . . . , v.sup.(M−1)].sup.T, following symbols, defined as an (L+1)×(L+K) filter matrix V.sub.L+1.sup.(l), are introduced as:
(41)
(42) wherein V.sub.L+1 is an M(L+1)×(L+K) matrix formed by stacking M filter matrices V.sub.L+1.sup.(l), (l=0, 1, . . . M−1), and defined as:
V.sub.L+1=[V.sub.L+1.sup.(0)T,V.sub.L+1.sup.(1)T, . . . V.sub.L+1.sup.(M−1)T].sup.T (12).
(43) Lemma 1: the following structural relationship holds:
v.sup.TH.sub.K=h.sup.TV.sub.L+1 (13).
(44) Obviously, the Lemma 1 holds because H.sub.K and V.sub.L+1 are row full rank and the number of columns of v.sup.T and that of h.sup.T are equal to the number of rows of H.sub.K and V.sub.L+1. According to the Lemma 1, a quadratic norm of a vector Ĝ.sub.i.sup.HH.sub.K is defined as:
∥Ĝ.sub.i.sup.HH.sub.K∥.sup.2=Ĝ.sub.i.sup.HH.sub.KĜ.sub.i=h.sup.Hĝ.sub.iĝ.sub.i.sup.Hh (14);
(45) where ĝ.sub.i is an M(L+1)×(L+K) filter matrix corresponding to a vector Ĝ.sub.i. The formula (10) may be represented as:
(46)
(47) According to the Theorem 1, if the quadratic form is constructed from an actual correlation matrix, then an actual channel impulse response coefficient is unique and therefore q(h)=0. Conversely, if only an estimated covariance matrix is available, a rank of the quadratic form is not equal to M(L+1)−1. Therefore, a nontrivial solution h can be solved out by minimizing the cost function q(h) under appropriate conditions. Typical constraints can be linear constraint and quadratic constraint.
(48) The quadratic constraint minimizes q(h) under a constraint of ∥h∥=1, to solve out a eigenvector corresponding to a minimum eigenvalue of the matrix Q as a solution. The quadratic form has a large calculation amount, but has a good stability.
(49) The linear constraint solves variable h of min[q(h)] under a constraint of c.sup.Hh=1, wherein c is an M(L+1)×1 vector. The solution at this point is proportional to Q.sup.−1c. The linear constraint has a relatively simple calculation, but has a poor stability due to a possibility of Q being singular.
(50) 4. Implementation of the Algorithm
(51) The channel impulse response coefficients can be calculated by quadratic minimization. The quadratic form as shown in the formula (15) may be equivalently represented in a form of a series of signal subspace eigenvectors, such as:
(52)
(53) where ŝ.sub.i represents an M(L+1)×(L+K) filter matrix corresponding to an eigenvector
(54)
The minimization of the formula (10) under the constraint of ∥h∥.sup.2=1 is equivalent to a maximization of the following formula (17):
{circumflex over ( )}q(h)=h.sup.H{circumflex over (Q)}h (17).
(55) The above formula (17) can be maximized by solving a maximum eigenvalue of {circumflex over (Q)}.
(56) The quadratic form with the noise as parameter has a same solution as the quadratic form with the signal as parameter under the constraint of unit norm, however calculation amount of the latter is less than the former, since the former needs to undergo MK−L−K times of summation operations and the latter only needs to undergo L+K times of summation operations.
(57) The theorem 1 points out that the channel impulse response coefficients can be determined by a noise subspace part of the correlation matrix R.sub.x. The following inference 1 can be derived according to the theorem 1.
(58) Inference 1: parameters of channels can be determined when only a few vectors of the noise subspace are known under appropriate conditions.
(59) The inference 1 is illustrated and proved below. It is proved that channel parameter vector H can be identified by utilizing a vector of the noise subspace for any K>L starting from a simplest case M=2 (i.e. two channels exist).
(60) The inference 1 can be proved according to the steps of: assuming K>L and G≠0 to be an arbitrary vector of the noise space, i.e. G=[G.sup.(0)T,G.sup.(1)T].sup.T for two channels M=2. Since G.sub.i.sup.HH.sub.K=O, 0≤i≤MK−L−K, i.e. the noise subspace is orthogonal to the column space of H.sub.K, it means:
G.sup.(0)(z)H.sup.(0)(z)+G.sup.(1)(z)H.sup.(1)(z)=0 (18);
(61) where G.sup.(i)(z), i=0, 1 is a z-transformation of a K×1 vector G.sup.(i), i=0, 1. H.sup.(0)(z) differs from H.sup.(1)(z) with only one constant in Formula (18) is proven by a contradiction analytical method as follows:
(62) It is assumed that H′.sub.K=[H′.sup.(0)T, H′.sup.(1)T is a 2(L+1)×1 vector, H′.sub.K represents a 2K×(L+K) filter matrix corresponding to H′. A vector G is orthogonal to column of H′.sub.K, if:
G.sup.(0)(z)H′.sup.(0)(z)+G.sup.(1)(z)H′.sup.(1)(z)=0 (19);
(63) where the condition that the formula (18) and the formula (19) cannot hold simultaneously is G.sup.(0)(z)H′.sup.(0)(z)+G.sup.(1)(z)H′.sup.(1)(z)=0. Since H.sub.K is assumed to be full-rank, H.sup.(0)(z) and H.sup.(1)(z) are relatively prime. With H.sup.(0)(z) being divided by H′.sup.(1)(z) and H.sup.(1)(z) being divided by H′.sup.(1)(z), and since the formula (18) and the formula (19) have the same dimensions, they must be proportional, that is:
(64)
(65) The inference 1 can be further proved in a manner similar to that described above for M≥2.
(66) The above inference 1 shows that the channel parameter vector H may be identified by using one of vectors in the noise subspaces, which means that H can be determined by using one item in the quadratic form.
(67) In fact, L+K linear equations including M(L+1) unknown parameters H will be generated if G.sub.0 is an MK×1 vector in the noise space of matrix R.sub.x, the vector G.sub.0 is orthogonal to the column space of H.sub.K, i.e. G.sub.0.sup.HH.sub.K=O or equivalently as H.sup.Hg.sub.0=0, where g.sub.0 is an M(L+1)×(L+K) filter matrix associated with G.sub.0. There will be no unique solution since the linear system cannot be determined if L+K<M (L+1). p such linear systems can be constructed by selecting p−1 additional vectors in noise subspaces. The system may be determine if the selected p meets p(L+K)≥M(L+1). The implementation and the process of the simulation experiment of the algorithm will be illustrated hereafter.
(68) (1) Implementation of the Algorithm
(69) Due to an overall average of the correlation matrix of communication signals being unknown in actual communication, it is replaced with a sample estimation. A specific implementation includes a step of:
(70) estimating the correlation matrix R.sub.x, wherein
(71)
and {tilde over (X)}.sub.i represents an observed sample value, thus R.sub.x can meet a Toeplitz structure;
(72) eigenvalue decomposition: decomposing the correlation matrix R.sub.x into a form of formula (8); and
(73) a calculation of a quadratic form: solving min[q(h)]=min[h.sup.rQh] according to formula (15), i.e., solving min(Q.sub.i)=min(ĝ.sub.iĝ.sub.i.sup.H), wherein h=Q.sup.−1c.
(74) Simulation Experiment:
(75) The algorithm is validated by using Monte Carlo simulation experiments. The transmitted signals adopt a series of independent 4QAM signals, and virtual channels adopt following channels:
(76)
(77) where the number of analog channels are M=4, a width of a time window K=10, an intersymbol interference length is L=4, an additional noise is Gaussian white noise, an output signal-to-noise ratio SNR=25 dB, and the numbers of signals in simulation are respectively 250, 500, 1000, 2000, and 5000. 100 times of simulations are performed by using Monte Carlo for each simulation to calculate an average deviation and a mean square error.
(78) (1) Average deviation: the average deviation is defined as:
(79)
(80) wherein H.sub.i
represents an average of the estimated values of H.sub.i under Monte Carlo simulation, and H.sub.i represents a i.sup.th component of an actual parameter vector.
(81) (2) Mean square error: being defined as an average of the mean square errors for each of all unknown parameters;
(82) (3) It can be seen from
(83) (4) The precision of the parameters estimated by using the noise vectors is related to the number of the noise vectors, it can be seen from
(84) (5) So the parameters of channels can be estimated effectively by using several noise vectors, under a condition of an appropriate number of samples, which fully demonstrates correctness and feasibility of the inference 1.
(85) The foregoing is only exemplary embodiments of the disclosure which is not intended to limit the disclosure in any way. Any simple amendments, equivalent variations and modifications of the embodiments according to the technical spirit of the disclosure are within the scope of the technical solution of the disclosure.