Adaptive Subchip Multipath Resolving for Wireless Location Systems

Reliable positioning of cellular users in a mobile environment requires accurate resolving of overlapping multipath components. However, this task is di ﬃ cult due to fast channel fading conditions and data ill-conditioning, which limit the performance of least-squares-based techniques. This paper develops two overlapping multipath resolving methods (adaptive and nonadaptive), and shows how the adaptive solution can be made robust to the above limitations by extracting and exploiting a priori information about the fading channel. Also the proposed techniques are extended when there are antenna arrays at the base station. Simulation results illustrate the performance of the proposed techniques.


INTRODUCTION
Wireless propagation suffers from multipath conditions.Under such conditions, the prompt ray may be succeeded by multipath components that arrive at the receiver within short delays.If this delay is smaller than the duration of the pulse shape used in the wireless system (i.e., the chip duration T c in CDMA systems), then the rays overlap.When this situation occurs, it results in significant errors in the estimation of the time and amplitude of arrival of the prompt ray. Figure 1 illustrates the combined impulse response of a two-ray channel using a conventional pulse shape in a CDMA IS-95 system in two situations.In the second situation, where the pulses overlap, the location of the peak is obviously delayed relative to the location of the prompt ray.Such errors in the time-of-arrival are particularly damaging in wireless location applications (a topic of significant relevance nowadays-see, e.g., [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17][18]).In these applications, small errors in the timeof-arrival can translate into many meters in terms of location inaccuracy.
There have been earlier studies in the literature on resolving overlapping multipath components (see, e.g., [19,20]).However, there are two sources of impairments that introduce significant errors into the resolution accuracy and which need special attention; these sources of error are particularly relevant in the context of mobile-positioning systems.The first impairment is due to the possibility of fast channel fading, which prohibits the use of long averaging intervals.This is because the estimation period in wireless location applications can reach up to a few seconds, which may cause the channel between the transmitter and the receiver to vary significantly during the estimation period, even for relatively slow channel variations.The second impairment is the possibility of noise enhancement, which occurs as a result of the ill-conditioning of the data matrices involved in most least-squares-based solutions.
In this paper, we develop an adaptive technique for resolving overlapping multipath components over fading channels for wireless location purposes.The technique is relatively robust to fast channel fading and data ill-conditioning.The following are the main contributions of this work.
(1) We first describe a framework for overlapped multipath resolving over fading channels via least squares.The framework will indicate why conventional least-squares techniques may fail for fading channels.
(2) We then point out the ill-conditioning problem that arises from using the pulse-shaping waveform deconvolution matrix.In order to avoid the possibility of noise enhancement as a result of this ill-conditioning, we show how to replace the least-squares operation by an adaptive solution.Still, while it avoids boosting up the noise, the adaptive filter solution might diverge or might converge slowly if not properly designed.To address this difficulty, we use a successive projection technique that incorporates into the design of the adaptive filter all available a priori channel information.
(3) We also describe a procedure for extracting a priori channel information and feeding it into the adaptive filter operation.(4) We then consider the case when there are multiple antennas at the base station and we extend the proposed algorithm for the single antenna case to systems with antenna arrays.

PROBLEM FORMULATION
Wireless propagation generally suffers from multipath conditions.A common model for the impulse response sequence of a multipath channel of length L is [21] where the {α l } and {x u l (n)} are, respectively, the unknown standard deviations (also referred to as gains) and the normalized fading amplitude coefficients (with unit variance); these coefficients model the Rayleigh fading effect of the channel.Several of the gains {α l } might be zero; and a nonzero gain at some l = l o would indicate the presence of a channel ray at the corresponding delay n = l o .Our strategy will be to estimate the gains {α l }, for all values of l, and then compare these values with a threshold.If any α l is lower than the threshold, then we set it to zero.In model (1), it is assumed that the sampling period for the sequence {h(n)} is a fraction of the chip duration, say for some integer N u > 1.In other words, the time variable

T s n
One-bit period (K chips) refers to multiples of T s and the superscript u in x u l denotes upsampling.By using an upsampled model for the channel impulse response, we will be able to resolve overlapping rays more accurately.Now consider the problem of estimating the gains {α l } from a received sequence {r(n)}, which is defined as follows:1 where is a known (upsampled) chipping sequence2 (its entries are 0 or ±1 when n is an integer multiple of N u ).The integer K denotes the processing gain of the communication system, that is, the ratio between the bit rate and the chip rate-see Figure 2.Moreover, {p(n)} P−1 n=0 is a known pulse-shape impulse response sequence, and v(n) is additive white Gaussian noise of variance σ 2 v .Let L r = P + KN u + L − 1 denote the total number of samples {r(n)}.To proceed with the analysis, we introduce the following assumption.
Assumption 1.The variations in the fading channel {x u l (n)} within the duration of the pulse-shaping waveform, p(n) (i.e., within a duration of P samples), are negligible.
This assumption is reasonable for cellular systems even for fast fading channels.For example, for an IS-95 pulse shaping waveform [22], the duration of the pulse shape is equal to 10T c , which corresponds to 8 microseconds.The autocorrelation function, R x u l (τ), of the fading sequence {x u l (n)}, at a time shift of 8 μs for a relatively fast mobile station (MS) moving at 60 mph and using a carrier frequency of f c = 900 MHz, is given from [21] by J o (2π80 × 8 × 10 −6 ) = 0.999994 ≈ 1.This high value for the autocorrelation between fading ray samples, {x u l (n)}, implies that they can be assumed to be constant within the assumed duration.Therefore, we may ignore variations in the coefficients {x u l (n)} within the pulse-shape duration.
Then the sequence that results from the first convolution p(n) [x u 0 (n)c u (n)] can be obtained as the matrix vector product: Likewise, the sequence that results from the second convolu- ] can be obtained as the matrix vector product and so on.We can represent the above results more compactly in matrix form as follows.Introduce the downsampled sequences: Then we obtain from (4) that where r is the received vector of length L r defined as r col r(0), r(1), . . ., r L r − 1 (11) and v is the noise vector defined as Moreover, C x is the KN u × L matrix defined as and h is the unknown path gain vector defined by In summary, the problem we are interested in is that of estimating h from the received vector r in (10), with the constraint that the matrix C x is not known completely since it depends on the unavailable quantities {x k ( j)}.To do so, we will exploit the statistical property of the fading quantities {x k ( j)}.

CONVENTIONAL MATCHED FILTERING
Let us examine first what happens if we correlate r(n) and c u (n) as in The result of this correlation (or despreading) operation is given (in vector form) by (1/K)C * r, where C is the following L r × L spreading code matrix: Then from (10), we get Nabil R. Yousef et al.

5
When K is large enough, and using the orthogonality property we obtain the approximation-see Appendix A: where A L is an L × L pulse-shaping convolution matrix similar to A, and X K is the L × L diagonal matrix x 0 ( j), . . ., Assuming ergodic processes, and taking the limit as K → ∞ of both sides of the above definition, we obtain lim Thus, unless the channel fading coefficients have static components, we get lim This result causes the output of the correlation process given by ( 17) to approach zero as K → ∞.Consequently, estimation techniques that are based on correlation (or matched filtering) will be unrobust when used to estimate the fading channels.This fact explains why it is difficult to obtain accurate location estimates using such techniques.

A PARTITIONED LEAST-SQUARES RECEIVER STRUCTURE
We now describe a technique for estimating h from (10) and which does not require knowledge of the {x k ( j)}.To begin with, note from (10) that if C x were known, then the leastsquares estimate for h could be found by solving which gives However, C x is not known since the {x k ( j)} themselves are not known.Thus we proceed instead as follows.
We first partition the received vector r in (10) into smaller vectors, say r m , of size NN u samples each (i.e., each r m contains N symbols with N u samples per symbol).Each r m will satisfy an equation of the form with {A m , C m x } similar to {A, C x } in (10) but of smaller dimensions, and where v m is defined by Then, in view of the earlier discussion, we are motivated to introduce the following algorithm.
(1) Partition the received vector r into M smaller vectors with NN u samples each, and such that the mth vector is given by Note that L r = MNN u .
(2) Introduce the NN u × L correlation (despreading) matrix and the L × L fading matrix X m : where n o = mN.Now N is usually small enough such that (m+1)N−1 j=no x l ( j) will not tend to zero and, hence, we will not be faced with the difficulty of having X m → 0, as was the case with X k (22).
(3) Multiply each vector r m from the left by (1/N)C * m , with m = 0, 1, . . ., M −1.The correlated (despreaded) output is denoted by At the same time N is large enough to get uncorrelated shifted spreading sequences, so that similar to (19), y m can be approximated by The despreaded vector y m can be used to estimate z m in the least-squares sense by solving (5) Introduce the vector β (averaged over all estimates z m ): For simplicity of notation, we will write |x| 2 to denote a vector whose individual entries are the squared norms of the entries of x: Using this notation, we can write The entries of β will be shown in the sequel to be related to estimates of the desired gains {α l }-see (49).

Parameter optimization and bias equalization
Assume that the length of the received data is large enough (L r → ∞).Then expression (36) becomes As M → ∞, the averaging process may be approximated by the expectation operation so that Using ( 30) and (31) gives which can rewritten as where and the pseudo-inverse matrix A † L is given by For mathematical tractability of the analysis, we introduce the following assumptions.
Assumption 2. The sequence {v (n)} is identically statistically independent (i.i.d) and is independent of each of the fading channel normalized gain sequences {x u l (n)}.
Although the sequence {v (n)} is not i.i.d, the assumption is a reasonable approximation in view of the fact that the entries of {v(n)} are i.i.d, and in view of the orthogonality of the spreading sequences.The argument in Appendix B, for example, shows that v (i) and v ( j) are uncorrelated for i = j.Assumption 2 is instead requiring the noises to be independent.It follows from (41 Assumption 3. The fading channel normalized amplitudes {x k ( j)} are statistically independent of each other.
This assumption is typical in the context of wireless channel modeling [21].Using (40), the elements of the vector β are individually given by Expanding, using Assumptions 2 and 3, and following the same procedure used in [23,24], it can be verified that

Doppler estimate
Optimal N Fading bias where B f (l) and B v (l) are, respectively, given by (45) In the above, R xk (i) is the autocorrelation function of each of the fading channel coefficients, that is, Expression (44) shows that β(l) includes a multiplicative fading bias B f (l) and an additive noise bias B v (l).Now consider the case of identical autocorrelation functions for all channel rays, say R x (i), and define the SNR gain This expression suggests an optimal choice for N by maximizing it with respect to N. A similar approach was used in [23,24] and N opt is found by solving the following equation: Once the {B f (l), B v (l)} have been estimated, they can be used to correct β(l) in order to estimate the channel gains { α l }: where The estimates B v (l) and B f (l) can be obtained by using the same procedure given in [23,24].Figure 3 shows the resulting multipath searcher.The Doppler estimate depicted in Figure 3 is required during the determination of R x (i) and, hence, N opt and the fading bias coefficient B f (l) [25].
Figure 4 shows the SNR gain for different values of Doppler frequencies.Moreover, N opt for different values of Doppler frequencies has been shown in Figure 5.

Difficulties
The main problem facing the least-squares multipath searcher of Figure 3 is the ill-conditioning of the pulse-shaping matrix A L , which increases with the sampling resolution.Figure 6 plots the condition number of the matrix A L (in dB) versus the oversampling factor N u .The ill-conditioning of A L results in noise enhancement, which in turn reduces the estimation accuracy.In the next sections, we explain how to use an adaptive filter solution in order to avoid the least-squares step and, more specifically, avoid the boosting up of the noise.In order to enhance the robustness of the adaptive solution, we will further show how to extract and incorporate into the design of the adaptive solution a priori knowledge about the multipath channel.

AN ADAPTIVE PROJECTION TECHNIQUE
We now describe an adaptive projection technique for channel estimation that exploits a priori information about the channel for enhanced accuracy.The technique replaces the least squares of Section 4 by an adaptive filter.The proposed method can be described as follows.
Recall that we need to solve least-squares problems of the form (32), that is, for successive values of m, where We will denote the entries of the successive y m by {d m (i)}.
Clearly, the solution of (51) can also be approximately attained by training an adaptive filter that uses the {d m (i)} as reference data and the rows of the L × L matrix A L as regression data.We will denote the rows of A L by {u i }.Since A L has only L rows, the adaptive filter is cycled repeatedly through these regression rows until sufficient convergence is obtained.In addition, it is explained in Appendix C how we can extract useful information about the channel such as its region of support (i.e., the region over which the channel taps are most likely to exist) and the largest amplitude that any of its peaks can achieve.This information can be exploited by the adaptive solution as explained below in order to enhance the accuracy and the resolution of the resulting multipath searcher.Thus the adaptive implementation can be described as follows.
(1) The received signal r(n) is applied to a bank of matched filters C * m in order to generate the vectors {y m }.  (2) A parallel-to-serial converter is applied to each y m in order to form the reference sequence {d m (i)}.
(3) An adaptive filter of weight vector w m i is used to estimate z m at the ith iteration (i.e., w m i is the estimate of z m at iteration i).The regression vector u i is obtained from the rows of A L .The adaptive filter is iterated repeatedly in a cyclic manner over the rows of A L until sufficient performance is attained.
(4) In addition, at every N p iterations, the weight vector of the adaptive filter is checked and, if necessary, a projection step P is performed in order to guarantee that the filter taps are consistent with the a priori information that is available about the channel taps.For instance, if we know that the channel has only two nonzero taps, then we zero out all taps except for the largest two taps (recall that since z m = X m h, and since X m is a diagonal matrix, then zero taps Feed information in the projection operation P Proceed as in the multipath searcher of Figure 3 α 0 α L−1 . . . in h translate into zero taps in the estimates of z m ).Specifically, the adaptive filter weight vector w m i is updated as follows: Here μ(i) is a step-size parameter, P refers to the projection procedure, and N p is an integer greater than or equal to one and less than or equal to the total number of iterations performed.
(5) The successive projections are based on information obtained from the upper branch of the block diagram in Figure 7.The first branch extracts information about the channel region of support and maximum amplitude.This information is extracted by noncoherently averaging the output of the matched filter bank to form J f (τ).The adaptive filter weight vector is successively projected onto the set of possible elements satisfying the constraints (e.g., tap locations and amplitudes should lie within the ranges specified by the a priori information).The adaptive filter weight vector is iterated till it reaches steady state.For instance, when the upper branch finds 3 taps, it gives a rough estimation for the location and amplitude of these taps.The projection scheme within the adaptive filter blocks checks the number of nonzero taps in w i , and forces the taps that are out of the detected range by the upper branch to zero.

Simulation results
The robustness of the proposed algorithm in resolving overlapping multipath components is tested by computer simulations.In the simulations, a typical IS-95 signal is generated, pulse shaped, and transmitted through various multipath channels.The total power gain of the channel components is normalized to unity.Figure 8 is a sample simulation that compares the output of the proposed adaptive algorithm to the output of the block least-squares multipath resolving technique of Section 4 for a two-ray fading multipath channel.The first plot shows the considered two-ray channel in the simulation.The second and third plots, respectively, show the output of block least-squares and block regularized least-squares stages.It is clear that both least-squares techniques lead to significant errors in the estimation of the time and amplitude of arrival of the first arriving ray.The last plot shows the output of the proposed estimation scheme.It is clear that the proposed algorithm is more accurate than leastsquares techniques.Here we may add that it was noted that the algorithm converges in around 30-50 runs.In this simulation, we have assumed 128 spreading sequences (K = 128), each chip is upsampled by order of 8 (N u = 8), the upsampled receiving vector is partitioned into 8 subblocks (M = 8), the receiving SNR before despreading is −15 dB and finally the adaptive filter step size is 0.7 (μ = 0.7).
Figure 9 shows the estimation time delay absolute error and amplitude mean square error of the prompt ray in overlapping multipath propagation scenarios versus the estimation period (T).The simulations are performed for various values of the maximum Doppler frequency ( f d ) and channel amplitude ratio.The results show a good ability of the proposed adaptive algorithm to resolve overlapping multipath components.In this simulation, we have assumed 128 spreading sequences (K = 128), each chip is upsampled by order of 8 (N u = 8), the upsampled receiving vector is partitioned into 8 subblocks (M = 8), the receiving SNR before despreading is −15 dB and the adaptive filter step size EURASIP Journal on Applied Signal  is 0.7 (μ = 0.7).Please note that different fading frequencies change the effective again after despreading according to (47).

RECEPTION WITH AN ANTENNA ARRAY
Using an antenna array at the base station can improve the location estimation by providing both the TOA and AOA information.An antenna array receiver integrates multiuser detection and beamforming with rake reception in order to  mitigate multiuser interference, cochannel interference and fading.
Thus consider an N a -element antenna array at the base station.In this case, the channel model ( 1) is replaced by where h(n) is now an N a × 1 vector, a(θ l ) is the N a × 1 array response as a function of the AOA of the lth multipath and it is given by a θ l = 1, e j2π(d/λ) sin(θl) , . . ., e j2π((M−1)d/λ) sin(θl) T .( Here, θ l is the AOA of the received signal over the lth multipath, d is the antenna spacing, and λ is the wavelength corresponding to the carrier frequency.Likewise, the received signal in (3) is replaced by where r(n) is now an N a × 1 vector.We can again use the arguments of Section 3 to replace (10) by where R is an L r × N a received matrix defined as R r 1 , r 2 , . . ., r Na (58) and r n is the received vector of length L r over the nth antenna array, that is, Moreover, V is the noise matrix where v n is the noise vector at the nth antenna array, and H is an LN a × N a Toeplitz path gain matrix whose first column is determined by Finally A θ is an L × LN a matrix that contains the array responses: where The problem we are interested in is that of estimating the {α l } from the received matrix R in (58).

The partitioned adaptive receiver
As in Section 4, we partition R into smaller matrices, R m , of size NN u × N a each.The matrix R m will then satisfy an equation of the form with {A m , C m x } similar to {A, C x } in (10) but of smaller dimensions, and where V m is defined by (66) Then, in view of the earlier discussion, we can use the same algorithm that we used in the case of single antenna.(68) The correlated (despreaded) output is denoted by When N is large enough, and similar to (31), y m can be approximated by The resulting signal y m in (70) is similar to the signal in (31), albeit with higher SNR due to the use of the antenna array.Therefore, the proposed estimation algorithm (32)-( 49) for the single antenna case can be used as well.
The system model for the resulting multiantenna adaptive receiver is illustrated in Figure 10.

Estimating the array response
We still need to estimate the array response matrix A θ .For the received signal R in (57) of size L r × N a , we define a correlation matrix, as in (17), as follows: where C was defined in (16) and K is the length of the spreading sequence.Now replace A θ H by Z, so that The least-square estimate of Z is given by Now, in order to estimate A θ from Z, we need an estimate of the channel matrix H.It can be estimated from (74) by noting that the matrix A θ,1 (the first L × L block of A θ ) is an identity matrix, so that Adaptive estimator given in Figure 4 α 0 . . .
Angle of arrival estimator a(θ 1 ) . . .a(θ L ) Array response matrix of each antenna where h is defined in (62) and h is an L × 1 vector that contains the first L elements of h.Moreover, v 1 and y 1 are the first column of V and Y, respectively.So h can now be estimated using (74) in the same manner as h was estimated from y m in (30) by using (49).Using h to create H, the leastsquares estimate of A θ can be obtained as (75)

Simulation results with antenna array
The robustness of the proposed algorithm in resolving overlapping multipath components when the base station has an array of antennas is tested by computer simulations.In the simulations, a typical IS-95 signal is generated, pulse shaped, and transmitted through various multipath channels.The total power gain of the channel components is normalized to unity.We have considered 4 antennas at the base station and Figure 11 compares the simulation results when there are multiple antennas and single antenna at the base station.In this simulation, we have assumed 128 spreading sequences (K = 128), each chip is upsampled by order of 8 (N u = 8), the upsampled receiving vector is partitioned into 8 subblocks (M = 8) and the adaptive filter step size is 0.7 (μ = 0.7).

CONCLUSIONS
This paper develops two overlapping multipath resolving methods (adaptive and nonadaptive), and illustrates how the adaptive solution can be made robust to fast channel fading and data ill-conditioning by extracting and exploiting a priori information about the channel.The proposed techniques are further extended to the case with antenna arrays at the base station.Simulation results illustrate the performance of the techniques.

APPENDICES A. PROOF OF (19)
To simplify (1/K)C * AC x , we start with the given A in (6) and express it as A Top(p), p col p(0), p(1), . . ., p(P − 1) , (A.1) where the notation Top(p) denotes the lower-triangular Toeplitz matrix determined by p.Let where Top(v) is n × l Toeplitz.Then (A.3) can be written as Due to the orthogonality property of the spreading sequences we have and, therefore, x l ( j) 0 . . .

C. EXTRACTING A PRIORI CHANNEL INFORMATION
In this appendix we explain how to extract useful a priori channel information from the received signal [26].This information is used in Section 5 by the adaptive searcher for resolving overlapping multipath components.
(1) A power delay profile (PDP) is evaluated as follows: (2) The region of support of the power delay profile, say R f , is determined by comparing the PDP with a threshold λ f .The region of support refers to the region of the delay (τ) that might contain significant multipath components: We restrict R f to the first continuous region of delays.
In other words, R f starts from the earliest delay that is higher than the threshold until the value of τ at which the PDP falls below the threshold.(3) The peak of the PDP is determined along with the delay that corresponds to the peak.Let τ f denote the delay of the peak of J f (τ): Moreover, let m f denote the value of the peak of J f (τ): The number of fading overlapping multipath components that exist in the region of support, R f , is determined by using the multipath detection algorithm of [26].Let the number of overlapping multipath components be denoted by O.
In summary, the following a priori information can be used in the multipath resolving stage.
(1) The delay of the ray to be resolved is confined to R f .
(2) The number of fading overlapping multipath components that exist in R f is equal to O. (3) The maximum amplitude of any ray in this region is less than or equal to the square root of the maximum value of J f (τ) after equalizing for the noise and fading biases that may arise in this value.This value is equal to C f (m f − B v ), where B v and C f are two noise and fading biases that can be calculated as described in (49)-(50).

10 f d = 20 f d = 30 f d = 40 f d = 50 f d = 80 Figure 4 :
Figure 4: SNR gain versus N for K = 256 and T c = 8.138 microseconds and for different Doppler frequencies.

Figure 7 :
Figure 7: An adaptive multipath resolving scheme using successive projections.
f d = 10 Hz f d = 40 Hz f d = 80 Hz

Figure 10 :
Figure 10: An adaptive multipath resolving scheme using an antenna array receiver.