An Adaptive Channel Estimation Algorithm Using Time-Frequency Polynomial Model for OFDM with Fading Multipath Channels

Orthogonal frequency division multiplexing (OFDM) is an e ﬀ ective technique for the future 3G communications because of its great immunity to impulse noise and intersymbol interference. The channel estimation is a crucial aspect in the design of OFDM systems. In this work, we propose a channel estimation algorithm based on a time-frequency polynomial model of the fading multipath channels. The algorithm exploits the correlation of the channel responses in both time and frequency domains and hence reduce more noise than the methods using only time or frequency polynomial model. The estimator is also more robust compared to the existing methods based on Fourier transform. The simulation shows that it has more than 5dB improvement in terms of mean-squared estimation error under some practical channel conditions. The algorithm needs little prior knowledge about the delay and fading properties of the channel. The algorithm can be implemented recursively and can adjust itself to follow the variation of the channel statistics.


INTRODUCTION
The 3G wireless communication system is the next generation mobile cellular system that aims to provide high rate data communications of a bit rate up to 2 Mbit/s. Among many technical challenges in this broadband system, the severe intersymbol interference (ISI) caused by multipath effect of wireless channels is an essential one. One effective technique to deal with this problem is the orthogonal frequency division multiplexing (OFDM) [1,2]. In OFDM systems, the entire bandwidth is partitioned into parallel subchannels by dividing the transmit data into several parallel low bit rate data streams to modulate the carriers corresponding to those subchannels. By doing so, the OFDM system has a relatively longer symbol duration, thus provides a great resistance to ISI and impulse noise. When the number of subchannels is large enough, the subchannels can be treated as independent of each other and only a one-tap equalizer is needed for each subchannel. Because of these advantages, OFDM has become a promising technique for broadband wireless communications.
Channel estimation is a key issue in a communication system, as is the case for the OFDM system. Without the knowledge of channel information, noncoherent detection, such as differential modulation, has to be used and results in some performance loss compared to the coherent detection. The channel estimation problem becomes more important for the 3G systems because many sophisticated signal processing techniques that require the knowledge of the channel information are expected to be used to meet the challenge of throughput and performance. For example, the independence of the subchannels in OFDM systems provides an easy way to optimize the transmitter design by adjusting the bit rate and transmit power across subchannels according to their channel conditions [3], which implies that the channel information has to be known at the transmitter.
The channel estimation problem is also more challenging in the 3G system because both the multipath effect and the fading effect have to be considered in this mobile broadband system. The important observation to solve the channel estimation problem in the OFDM systems is that the fading multipath channel in the 3G system is correlated in both time and frequency domain, even though subchannels are treated independently when performing the signal detection. The channel estimation algorithms should exploit such correlation to improve the accuracy of the estimation. Van de Beek et al. [4] tried to exploit the correlation of the channel parameters in frequency domain while Mignone and Morello [5] used the correlation in time domain. Li et al. [6] considered the correlation in both time and frequency domains. The estimators designed in these literatures are all Fouriertransform-based approaches, which implicitly assumed that the channel power spectrum can be viewed as band limited. The assumption is true when we consider the ensemble statistics. However, in practice, we can only get finite discrete samples of the channel response of the time varying channel. The leakage can be very severe and then degrade the performance dramatically.
In this work, we consider the problem from another point of view. Because of the correlation of the fading multipath channel, it can be viewed as a smoothly varying function of both time and frequency. It has been stated in the approximation theory that such a smoothly varying function can be approximated by a set of basis functions [7], for example, the polynomial basis [8]. Borah and Hart [9,10] used the time domain polynomial approximation while Luise et al. in [11] used the frequency domain polynomial approximation. However, the channel responses used in coherent detection of OFDM are located in the time-frequency plane. Therefore, it is naturally to exploit the channel correlation in both time and frequency domains using a time-frequency polynomial model. The noise can then be greatly suppressed by estimating a smaller number of coefficients of the basis functions over a large number of observations. Moreover, it also make the estimator design more flexible and robust to the variation of channel statistics. We can also view Fourier transform as a type of model basis function and hence Fourier-transformbased method is the same type of method as the polynomialmodel-based channel estimation scheme but with different model accuracy and different noise reduction capability. These two methods compared, the model error of the Fourier basis is very sensitive to the channel statistics and works only for some very specific system parameters and channel statistics. On the contrary, the polynomial-model-based method performs more consistently and robustly for variety of channels.
A key problem in using the polynomial model to estimate the channel responses is to decide the model order and timefrequency window dimensions of observations. The model approximation error of polynomial model decreases when increasing model order or decreasing the window dimensions. On the other hand, the noise is reduced more when decreasing the model order or increasing the window dimensions. It is important to reach a tradeoff between the model error and noise reduction. In this paper, we propose an adaptive algorithm that adjusts the window dimensions to balance the tradeoff. With this adaptive algorithm, the channel correlation function or the fading and delay characteristics are no longer that essential in the design of the channel estimator. The estimator can adapt its settings to the variation of the channel statistics.
. . The rest of the paper is organized as follows. First, we introduce the OFDM system in Section 2 and the fading multipath channel in Section 3.1. Then, we discuss the time-frequency polynomial model Section 3.2 and derive the corresponding recursive channel estimation algorithm in Section 4. The performance analysis is discussed on the general-model-based estimation approach in Section 5. Then the window dimension adaptive algorithm is derived in Section 6 based on the performance analysis. Finally, the simulation results are presented to demonstrate the performance in Section 7 and the conclusion is drawn in Section 8. Figures 1a and 1b show the transmitter and receiver of an OFDM system, respectively. The OFDM system divides the whole bandwidth B d into m subchannels by buffering the input data to blocks, and then partitions the block into m lower rate bit streams. In most of OFDM systems, the subchannels are divided evenly, the bandwidth of the subchannels or the rate of the bit streams is ∆ f = B d /m. The bit streams may contain different amount of bits and use different transmit energy according to the channel condition. The bit and energy allocation is done by a loading algorithm. Then the bit streams are mapped to some complex constellation points X i,k , i = 0, . . . , m−1 at the kth block. The modulation is then implemented by m-point inverse discrete Fourier transform (IDFT). Then the modulated data go through P/S converter to form the serial data x i,k . A cyclic prefix which is constructed using the last v samples of x i,k 's is inserted before sending the x i,k 's to the channel. Now it follows that the symbol duration is m/B d , however, the actual block duration is  T f = (m + v)/B d . For a system with B d = 800 kHz, m = 128, and v = 16, the block duration is T f = 180 microseconds. Such a system is used in the rest of this paper.

OFDM SYSTEMS
At the receiver, the prefix part is discarded. The demodulation is performed by the discrete Fourier transform (DFT) operation. If the cyclic prefix is long enough, then the interference between two OFDM blocks is eliminated and the subchannels can be viewed as independent of each other, that is, the demodulated data Y i,k can be expressed as where H i,k is the channel frequency response at i∆ f of kth block and N i,k is the corresponding channel noise that is assumed to be white Gaussian process with zero mean and variance σ 2 . Because of the simple relation of (1), only a one-tap equalizer is needed for each subchannel at the receiver, that is,X where the equalizer coefficient W i,k is some function of H i,k . For example, the zero-forcing equalizer is constructed as Then the decision or decoding is made upon X i,k .

Fading multipath channel
In a mobile broadband wireless communication system such as 3G, the transmission is impaired by both fading that is due to the mobility, and multipath that is due to the wide bandwidth. This fading multipath channel has long been known to be modeled as a time-varying linear filter [12], where γ i (t)'s are independent complex Gaussian processes with zero mean and variance p i 's. For OFDM systems, we can assume that the channel is time varying for different blocks but time-invariant within one block. The channel frequency response H i,k 's are samples of the continuous channel response The correlation function of For the Rayleigh fading channel [12], r t (t) = J 0 (2π f D t) and r f ( f ) = i p i e − j2π∆ f τi with J 0 (·) denoting the zero-order Bessel function, f D being the Doppler frequency describing the channel variation along t, and p i 's together with τ i 's being delay profiles describing the channel dispersion which is also often characterized by the maximum delay spread T d = max i τ i . Three types of delay profiles are used in this work, TU, HT, and 2-ray [13]. The TU and HT delay profiles are shown in Figure 2. The 2-ray delay profile has two equal power paths and the delay between two paths is T d . We also assume that the channel is normalized in our simulation, that is, i p i = 1. From (5), the power spectrum of the channel response is Because of the physical mechanism of the propagation, the channel varies smoothly and most of the energy is concen-trated in a finite bandwidth in both time and frequency domains. The bandwidths are f D for S t (ξ) and T d for S f (ν), respectively.

Time-frequency polynomial channel model of OFDM systems
We know from the approximation theory [7,8] that the smoothly varying channel responses can be approximated by projecting to a finite set of basis functions. In [14], it was shown that the channel responses in a small time domain window around a center point k 0 of dimension 2K + 1 can be closely approximated by a small set of polynomial basis functions, that is, For such an approximation, it can be proved that the mean-squared model error is bounded by Here we assume that fD 0 S t (ξ) dξ = 1. It can be seen from (10) that the sufficient condition for this error to converge to zero is Similarly, if the channel delay spread T d satisfies T d ∆ f 1, which means that the frequency variation of H i,k 's is smooth enough along frequencies, then can be approximated by the polynomial bases, that is, where with The mean-squared model error of this approximation is bounded by The time domain expansion (7) is used for channel estimation in [10], while the frequency domain expansion (10) is adopted in [11]. For a given channel, the selection of the above two types of expansions depends on the channel statistics, f D and T d , and the system parameters, T f and ∆ f . Moreover, it is naturally to expand the channel responses in both time and frequency domains [15] for the OFDM system, since its signal is distributed in a time frequency plane. The expansion in the time-frequency window of dimensions (2I + 1) × (2K + 1) around i 0 and k 0 is with The mean-squared model error is then bounded by Without loss of generality, assuming M = N and using the multipath Rayleigh fading channel described in Section 3.1, we can show that [15] Again the sufficient conditions for convergence of the above expansion are f D T f 1 and T d ∆ f 1. It is noticed that these two conditions are usually satisfied in a practical OFDM systems. When f D T f > 1, the block duration is so long or the channel changes are so fast that the channel cannot be viewed as invariant during one block and the system suffers large interchannel interference (ICI). On the other hand, when T d ∆ f > 1, the block duration is so short or the channel dispersion is so large that the subchannels can no longer be treated independently and the system would suffer both ISI and ICI. The OFDM system cannot work in either case. Hence, it is reasonable to assume that both conditions are satisfied in a well designed OFDM system. Now we take a close look at the upper bound of the model error. Suppose that the length of the cyclic prefix v can be ignored compared to the number of the subchannels. Then the first term in (16) is determined by The third term is actually determined by f D T d and is much smaller than the first two terms, since they both are smaller than one. For the first two terms, when m is large, the first term is dominating, then we should choose smaller K or larger M. While m is small, then the second term is dominating and I should be smaller or the model order N should be larger.
If the Doppler frequency f D , maximum delay T d , and bandwidth B d are fixed, we can adjust the window dimensions according to m to keep the time-frequency model error to certain level but we still have a small MN/IK. However, if only time or frequency domain expansion is used, the model error cannot be adjusted to maintain a small level with the same M/K or N/I when the number of subchannels m varies. Figure 3 shows the upper bound of mean the squared model error with f D T f = T d ∆ f = 10 −2 and I = K = 5 according to model order M. It shows that the model error is under −40 dB as the model order is 3. This means that we only need to estimate 9 model coefficients to get the 121 channel responses. In this figure, we also show the residual noise for SNR of 10 dB and 20 dB. It shows that the noise can be greatly reduced with very small penalty on model error. Moreover, such a model approximation does not need to know the actual channel correlation function.

Estimator structure
The channel estimation problem in OFDM systems is to estimate the channel response H i,k based on the transmitted signal X i,k and the received signal Y i,k . The information of the transmitted signal X i,k 's is obtained either from training or from detected feedback. In OFDM systems, an instantaneous estimate can be easily constructed Then suppose that we have chosen the model order and window dimensions such that the model error is small and can be ignored, we can approximate (13) in a matrix form where q I,K 00 · · · q I,K 0,N−1 · · · q I,K M−1,0 · · · q I,K Definẽ where Using least square (LS) methods, we can get the estimation of the coefficients of the polynomial basis from the instantaneous estimatesb where Q † M,N (I, K) is the pseudoinverse of Q M,N (I, K). The channel estimation then can be constructed aŝ where Usually, we fix the value of i − i 0 and k − k 0 , that is, we fix the point of estimation inside the window and slide the window to get all the estimations. Then the estimator can be viewed as a two-dimensional filtering process. Arranging the instantaneous estimation inside the window into a matrix form, Then the estimation iŝ where q N (i) = i 0 i 1 · · · i N−1 T , The estimator structure is shown in Figure 4. The coefficients of the frequency domain filter are Q † N (I)q N (i − i 0 ) and the coefficients of the time domain filter are Q † M (K)q M (k − k 0 ).

Recursive algorithms
The two-dimensional filter can actually be implemented recursively in time and frequency domains, respectively. Define the basis functions Q N (I) as Q N (I)= q N (−I + 1) · · · q N (0) · · · q N (I) q N (I + 1) T .

(28)
The basis Q N (I) and Q N (I) are actually homomorphic to each other, that is, there is an invertible matrix R such that This means that instead of using Q N (I) as basis, we can use Q N (I) as the basis to construct the estimator, that is,H Substitute (30) into (26), we can estimate the channel using Q † n (I). Then the core of the recursive algorithm is to calculate Q † N (I) from Q † N (I) iteratively. Let P f = (Q T N (I)Q N (I)) −1 and P t = (Q T M (K)Q M (K)) −1 . At initialization, we estimate model coefficientsb f (k) or b t (i) regarding to P f or P t over a window of dimension 2I +1 or 2K + 1. Similar to the recursive least square (RLS) algorithm, using the matrix inverse lemma [16], we can calculate P + f = (Q T N (I + 1)Q N (I + 1)) −1 or P + t = (Q T M (K + 1)Q M (K + 1)) −1 and then the corresponding model coefficientsb + f or b + t recursively over the window of dimension 2I + 2 or 2K + 2 during the updating process. After that apply the matrix inverse lemma again, we can calculate P − f = (Q T N (I)Q N (I)) −1 or P − t = (Q T M (K)Q M (K)) −1 and the corresponding model coefficientsb − f orb − t from P + f or P + t over the window of dimension 2I + 1 or 2K + 1 during the downdating process. Then according to (30), the channel can be estimated aŝ As this recursive process going on, the basis function becomes q T M (k + l − k 0 ) where l is the index of the iteration. The dynamic range of such a basis function may become so large that it will affect the numerical stability of the algorithm. Therefore, regularization using R f or R t should be used periodically to scale the basis back to Q N (I) or Q M (K). The frequency domain and time domain recursive algorithms are summarized in Algorithms 1 and 2, respectively. The matrices K + f and K − f or K + t and K − t are the corresponding gain matrices in updating and downdating. The two-dimensional filter in the tables is implemented first by frequency domain filtering then by time domain filtering. The order can be switched. In that way, the input in Algorithm 2 are instantaneous estimates while the inputs in Algorithm 1 are the outputs of the time domain filters of Algorithm 2. It is also noted that the order of downdating and updating can be switched, too.
In both tables, Updating: with the new inputH I+i0+1,k , calculate , : , The recursive algorithm needs less calculation compared to direct computation of the product of pseudoinverse when the window dimensions are much larger than the model order. Many fast algorithms of recursive least square (RLS) can be used for the practical implementation of such a recursive algorithm [16]. It also provides an easy way to adjust the window dimensions for the implementation of the adaptive algorithm.

PERFORMANCE ANALYSIS
Suppose that the channel can be modeled by some basis function, that is, a set of channel responses H can be projected to a set of basis function Q and the coefficients of the basis functions are b, that is, and the other is related to the residual noise, that is, Since R H is a Toeplitz matrix, it can be decomposed as where Λ is a diagonal matrix with eigenvalues of R H on its diagonal. If Q = U 1 , then the model error is zero. This means that the optimal function basis, which we can find in terms of model accuracy, is the eigenbasis U 1 . However, it requires the knowledge about the statistics of the channel responses.
In some special cases, we can easily find some specific function bases that can diagonalize R H without actually knowing R H . For example, if H is the channel response for one OFDM block with all the delay paths, τ i 's, at the sampling instance of the OFDM system, then such an optimal function basis is the DFT matrix [6]. However, in most of practical situations, the channel delay profiles do not satisfy this condition. Therefore, using DFT matrix may cause severe leakage problem and incur a large model error. The average energy of the residual noise over the entire estimation window can be calculated as follows: The average mean-squared error over the whole estimation window is actually lower bounded by (41). The lower bound is achieved when Q = U 1 .
Although the average energy of the residual noise maintains the same once the data length and model order is fixed, the estimation error inside the window is often distributed unevenly and differently for different basis functions. For the polynomial model, the estimation error is the least at the center point of the window and larger at the edge. Therefore, we prefer to choose the center of the window to get a better performance. However, along the time domain, we can only choose the end point to get a causal filter.

OPTIMAL MODEL PARAMETERS ADAPTATION
With estimation point chosen at the center of the frequency domain window and end point at the time domain window, the estimation error from (23) becomes where the model error is   The residual noise is reduced more when the model order M × N becomes small or the window dimension I × K becomes large. However, the model error will increase in this case. With fixed polynomial model order M and N, the optimal window dimension is obtained by min I,K Usually, there are several local minima in this optimization problem. Considering the computational complexity, we would prefer the one with small I × K.
In order to adaptively adjust the window dimensions we need to know the estimation error. Since the actual channel responses are not known, we have to estimate the estimation error using the instantaneous estimates and the final estimates. Suppose that the noise statistics is known, we can calculate the estimated estimation error aŝ Using this approximation, the window dimension adaptive algorithm for the optimization of (45) is given in Algorithm 3.
If the recursive algorithm in Algorithms 1 and 2 is adopted, the window adaptation can be implemented easily. We just eliminate one downdating when increasing the window dimension, or eliminate one updating when decreasing window dimension.
One important problem in the adaptive algorithm is to determine the threshold f th and t th . With large threshold, the algorithm converges faster, but with larger deviation. Especially when the local minima are located closely, the large threshold may result in unstable convergence. Hence, it would be preferred to use smaller thresholds here.

SIMULATION RESULTS
The OFDM system used in the simulations is the system introduced in Section 2. QPSK modulation is used throughout all subchannels. Figure 5 shows the mean-squared estimation error and the symbol error rate (SER) comparison of the al- is better than that of using only frequency domain expansion or using only time domain expansion in both cases. However, in the first case, the delay spread is smaller while the Doppler is larger, then the channel responses have more correlation in the frequency domain than in the time domain. Therefore, we use larger frequency domain window to exploit the frequency domain correlation. In the second case, the delay spread is larger while the Doppler frequency is smaller, then the channel responses have more correlation in the time domain and we use a larger time domain window to exploit it. It is shown that we have to use different time and frequency estimator to best exploit the channel correlations for different channels. Using only time or frequency domain scheme is not enough. Figure 6 shows the estimation error under different delay profiles with Doppler frequency of 40 Hz. Figure 6a shows the estimation error with TU delay profile and 2-ray delay profile of T d = 5 microseconds, which is the maximal delay spread of TU while Figure 6b shows the estimation error with HT delay profile and 2-ray delay profile of T d = 17.2 microseconds which is the maximal delay spread of HT. We also compared the results using the Fourier-transform-based method of [6]. We can see that for TU or HT, the proposed algorithm performs much better than the Fourier-transform-based method. However, for 2-ray channel with T d = 5 microseconds, the Fouriertransform-based method performs the best. The reason is that T d = 5 microseconds is an integer multiplication of the sampling period of the OFDM system, which is t s = 1/800 KHz = 1.25 microseconds. The impulse response of this 2-ray channel has energy only at the sampling instance  of the OFDM system, hence there is no leakage or model error using Fourier transform, which is used as frequency domain estimator in [6]. In this case, the Fourier-transformbased method actually provides a minimum mean-squared error estimator. Unfortunately, in the practice, such a case is quite unlikely especially for the time-varying channel. It is show that for TU or HT delay profiles and 2-ray with T d = 17.2 microseconds, there is great amount of leakage  using Fourier-transform-based on the sampling frequency of the OFDM system. The leakage greatly degrades the performance of the Fourier-transform-based method. In contrast, the polynomial-model-based method performs consistently for the channels with same maximal delay spread and hence is more robust to the channel statistics. This is because the model errors are bounded by the same bound for the channels with the same T d and f D as stated in Section 3.2. Therefore, the performance of the polynomial-model-based channel estimation is not sensitive to the specific correlation functions of the channels with the same Doppler frequency and maximum delay spread. Figure 7 shows the mean-squared estimation error at SNR of 20 dB with different delay spread of a 2-ray channel. It further demonstrates the robustness of the polynomialmodel-based method compared to the Fourier-transformbased method. The Fourier-transform-based method performs better only when the delay spread is at the sampling instance of the system. For most of the cases, it performs poorly. However, for the polynomial-model-based method, it performs consistently and outperforms the Fouriertransform-based method most of the time. Figure 8 shows the window dimension adaptation. The window dimension variation is shown in Figures 8a and 8b. The estimation error is shown in Figures 8c and 8d. Two cases with different initial conditions are simulated, which are shown in Figures 8a and 8c and Figures 8b and 8d, respectively. In Figures 8a and 8c, the window dimension is 5 × 5 at the beginning, while in Figures 8b and 8d, it is 20 × 20. In both cases, after about 100 iterations, the algorithm converges to a window dimension of 12 × 10 and an estimation error under −26 dB. However, as mentioned in Section 6, smaller window dimensions are preferred for the sake of the computation complexity. With this adaptation algorithm, the polynomial-model-based method is not only robust to the specific correlation of the channel variation and dispersion, but also robust to T d and f D and can follow the variation of the statistics of the channel. Moreover, in the previous simulation, fixed window dimensions are used, by applying this window dimension adaptation algorithm, the performance in Figure 6 can be further improved.

CONCLUSIONS
In this work, we proposed a channel estimation algorithm for the OFDM system with fading multipath channels, which is suitable for the applications in 3G wireless communications. The algorithm is based on the time-frequency polynomial model that exploits the correlation of the channel responses in both time and frequency domains. The channel response is approximated by a small number of time-frequency polynomial basis functions and estimated by first estimating the coefficients of the bases. The residual noise is significantly reduced in this way, compared to the results when approximation is only done either in time or frequency domain, and the estimator design is more flexible. Therefore, the approach is more robust to the channel statistics and system parameters than the existing Fourier-transform-based method. It does not require the delay profiles to be integer multiples of the system sampling period. Moreover, the algorithm can be implemented recursively and can adjust the model parameters adaptively to the delay and fading characteristics.