Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Channel estimation techniques for OFDM systems in dispersive time-varying channels Song, Xuegui 2009

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata

Download

Media
[if-you-see-this-DO-NOT-CLICK]
[if-you-see-this-DO-NOT-CLICK]
ubc_2009_fall_Song_Xuegui.pdf [ 615.06kB ]
[if-you-see-this-DO-NOT-CLICK]
Metadata
JSON: 1.0067250.json
JSON-LD: 1.0067250+ld.json
RDF/XML (Pretty): 1.0067250.xml
RDF/JSON: 1.0067250+rdf.json
Turtle: 1.0067250+rdf-turtle.txt
N-Triples: 1.0067250+rdf-ntriples.txt
Original Record: 1.0067250 +original-record.json
Full Text
1.0067250.txt
Citation
1.0067250.ris

Full Text

C HANNEL E STIMATION T ECHNIQUES FOR OFDM S YSTEMS IN D ISPERSIVE T IME -VARYING C HANNELS by  X UEGUI S ONG B.Eng., Northwestern Polytechnical University, P. R. China, 2003  A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF  MASTER OF APPLIED SCIENCE in The College of Graduate Studies (Electrical and Computer Engineering)  THE UNIVERSITY OF BRITISH COLUMBIA (O KANAGAN) May 2009 c Xuegui Song, 2009  Abstract Coherent modulation is more effective than differential modulation for orthogonal frequency division multiplexing (OFDM) systems requiring high data rate and spectral efficiency. Channel estimation is therefore an integral part of the receiver design. In this thesis, two iterative maximum likelihood based channel estimation algorithms are proposed for an OFDM system in dispersive time-varying channels. A multipath channel model is proposed for OFDM uplink transmission in macrocellular systems. The multipath fading channel is modeled such that the channel state can be determined by estimating the unknown channel parameters. A second-order Taylor series expansion is adopted to simplify the channel estimation problem. Based on the system model, an iterative maximum likelihood based algorithm is first proposed to estimate the discrete-time channel parameters. The mean square error performance of the proposed algorithm is analyzed using a small perturbation technique. Based on a convergence rate analysis, an improved iterative maximum likelihood based channel estimation algorithm is presented using a successive overrelaxation method. Numerical experiments are performed to confirm the theoretical analyses and show the improvement in convergence rate of the improved algorithm.  ii  Table of Contents Abstract  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  ii  Table of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  iii  List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  v  List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  vi  List of Acronyms  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viii  Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  x  1  Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  1  1.1  Background and Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  1  1.2  Literature Review . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  2  1.3  Thesis Outline and Contributions  . . . . . . . . . . . . . . . . . . . . . . . . . .  4  Wireless Channels and OFDM Technique . . . . . . . . . . . . . . . . . . . . . . . .  6  2.1  6  2  2.2  Wireless Channels  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  2.1.1  Large-Scale Fading  . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  6  2.1.2  Small-Scale Fading  . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  8  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  9  OFDM Technique 2.2.1  Basic Concept . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10  2.2.2  OFDM Implementation Using FFT . . . . . . . . . . . . . . . . . . . . . 11  2.2.3  Advantages and Drawbacks of OFDM  . . . . . . . . . . . . . . . . . . . 13  iii  Table of Contents 3  ML-based Channel Estimation for OFDM Systems in Dispersive Time-Varying Channels . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15  4  5  3.1  System Model  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15  3.2  Formulation of Channel Estimation Problem  3.3  An Iterative Channel Estimation Algorithm . . . . . . . . . . . . . . . . . . . . . 22  3.4  Performance Analysis of the Channel Estimation Algorithm . . . . . . . . . . . . 23  3.5  Numerical Results and Discussions  . . . . . . . . . . . . . . . . . . . . 16  . . . . . . . . . . . . . . . . . . . . . . . . . 28  Convergence Study of ML-based Channel Estimation Algorithms . . . . . . . . . . 35 4.1  Convergence Performance Analysis of the Proposed Algorithm  4.2  An Improved Channel Estimation Algorithm  4.3  Numerical Results and Discussions  Conclusions  . . . . . . . . . . 35  . . . . . . . . . . . . . . . . . . . . 41  . . . . . . . . . . . . . . . . . . . . . . . . . 41  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50  5.1  Summary of Contributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50  5.2  Future work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51  Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52  Appendices A Derivation of the Discrete CIR Model . . . . . . . . . . . . . . . . . . . . . . . . . . 56 B Derivation of (3.20) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58 C Derivation of (3.30c) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61 D Derivation of the Absolute Value for the Largest Eigenvalue of Mη−1 Nη  . . . . . . . 63  E Matlab Code . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65  iv  List of Tables 3.1  Algorithm 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22  4.1  Algorithm 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41  v  List of Figures 2.1  Signal power fluctuation versus range in wireless channels. . . . . . . . . . . . . .  7  2.2  Basic principle of OFDM modulation. . . . . . . . . . . . . . . . . . . . . . . . . 12  2.3  Basic principle of OFDM demodulation. . . . . . . . . . . . . . . . . . . . . . . . 12  2.4  OFDM modulation by means of IFFT processing. . . . . . . . . . . . . . . . . . . 13  2.5  OFDM demodulation by means of FFT processing. . . . . . . . . . . . . . . . . . 14  3.1  An MSE comparison for three different approximation orders. . . . . . . . . . . . 20  3.2  An MSE performance comparison between Algorithm 1 and Moose’s estimator for fˆl when fl = 0.05, δh = 1% and δ f = 1%. . . . . . . . . . . . . . . . . . . . . . . 30  3.3  An average MSE performance comparison between the simulated results and the theoretical values for hˆ l of Algorithm 1 when δh = 1% and δ f = 1%. . . . . . . . . 31  3.4  An MSE performance comparison between the simulated results and the theoretical values for fˆl of Algorithm 1 when fl = 0.02, δh = 1% and δ f = 1%. . . . . . . 32  3.5  An MSE performance comparison between the simulated results and the theoretical values for fˆl of Algorithm 1 when fl = 0.04, δh = 1% and δ f = 1%. . . . . . . 33  3.6  An MSE performance comparison between the simulated results and the theoretical values for fˆl of Algorithm 1 when fl = 0.06, δh = 1% and δ f = 1%. . . . . . . 34  4.1  An average MSE performance comparison between Algorithm 1 and Algorithm 2, along with the theoretical values for hˆ l , when δh = 1% and δ f = 1%. . . . . . . . . 43  4.2  An MSE performance comparison between Algorithm 1 and Algorithm 2, along with the theoretical values for fˆl , when fl = 0.02, δh = 1% and δ f = 1%. . . . . . 44  vi  List of Figures 4.3  An MSE performance comparison between Algorithm 1 and Algorithm 2, along with the theoretical values for fˆl , when fl = 0.04, δh = 1% and δ f = 1%. . . . . . 45  4.4  An MSE performance comparison between Algorithm 1 and Algorithm 2, along with the theoretical values for fˆl , when fl = 0.06, δh = 1% and δ f = 1%. . . . . . 46  4.5  An average number of iterations comparison between Algorithm 1 and Algorithm 2 when η = 4/3, δh = 1% and δ f = 1%. . . . . . . . . . . . . . . . . . . . . . . . 47  4.6  The average MSE performance versus number of iterations for hˆ l of Algorithm 1 and Algorithm 2 when η = 4/3 and SNR = 25. . . . . . . . . . . . . . . . . . . . 48  4.7  The average MSE performance versus number of iterations for fˆl of Algorithm 1 and Algorithm 2 when η = 4/3 and SNR = 25 dB. . . . . . . . . . . . . . . . . . 49  vii  List of Acronyms Acronyms  Definitions  3GPP  Third Generation Partnership Project  BER  Bit-Error Rate  BPSK  Binary Phase-Shift Keying  CIR  Channel Impulse Response  CRB  Cramer-Rao Bound  DAB  Digital Audio Broadcasting  DFT  Discrete Fourier Transform  DVB  Digital Video Broadcasting  FFT  Fast Fourier Transform  ICI  Inter-Carrier Interference  IDFT  Inverse Discrete Fourier Transform  IFFT  Inverse Fast Fourier Transform  ISI  Intersymbol Interference  LOS  Line-of-Sight  LS  Least Square  LTE  Long Term Evolution  ML  Maximum Likelihood  MMSE  Minimum Mean Square Error  MSE  Mean Square Error  OFDM  Orthogonal Frequency Division Multiplexing  P/S  Parallel-to-Serial  PAPR  Peak-to-Average Power Ratio viii  List of Acronyms RF  Radio Frequency  S/P  Serial-to-Parallel  SC-FDMA  Single-Carrier Frequency Division Multiple Access  SNR  Signal-to-Noise Ratio  SOR  Successive Overrelaxation  WLAN  Wireless Local Area Network  WMAN  Wireless Metropolitan Access Network  ix  Acknowledgments I am deeply grateful to my supervisor Dr. Julian Cheng for his enthusiasm, guidance, advice, encouragement, and support. I will continue to be influenced by his rigorous scholarship, clarity in thinking, and professional integrity. I would like to thank Dr. Zheng Du from Huawei Technology Co. and Dr. Norman C. Beaulieu from University of Alberta, Edmonton for their feedbacks, constructive comments, and valuable suggestions on my research work. I would like to express my thanks to Dr. Paramjit Gill for his willingness to serve as my external examiner. I would also like to thank Dr. Richard Klukas and Dr. Jonathan Holzman for serving on the committee. I really appreciate their valuable time and constructive comments on my thesis. I owe many people for their generosity and support during my two-year study at the University of British Columbia. I would like to thank my dear colleagues Chiun-Shen Liao, Junfeng Zhao, James Jianchen Sun, Mingbo Niu, Nick Kuan-Hsiang Huang, Ning Wang, Xian Jin, and Yeyuan Xiao for sharing their academic experiences and constructive viewpoints generously with me during our discussions. I would also like to thank my dear friends Eric Jinglong Ma, Junning Li, Pan Luo, Yi Zhang, and Yu Hui for sharing in my excitement and encouraging me when I was in frustration during this journey. I would also like to thank all my friends in Vancouver, Beijing, Chengdu, and Xi’an. Finally, I would like to thank my parents for their patience, understanding and support over all these years. All my achievements would not have been possible without their constant encouragement and support.  x  Chapter 1 Introduction 1.1 Background and Motivation The earliest prototype of wireless communications can be traced back to the pre-industrial age. In those years, information was transmitted over line-of-sight (LOS) distances using smoke signals, torch signaling, flashing mirrors, signal flares, or semaphore flags. However, true wireless communications, where information is transmitted in terms of electromagnetic waves through complex physical mediums, started from the first experiments with radio communication by M. G. Marconi in the 1890s. Since then, wireless communication systems have been developing and evolving with a rapid pace. In the intervening hundred years, many types of wireless communication systems have flourished, and often later disappeared. By far the most successful wireless communication system has been the cellular system. It is reported that the number of mobile phone subscribers is on the order of 4 billion currently and will rise to 5.6 billion in 2013 worldwide. Undoubtedly, cellular phones have radically changed the life of people and become a critical business tool in most countries. In the last two decades, with the explosive growth of wireless communications, wireless services have migrated from the conventional voice-centric services to data-centric services. Therefore, one main target for modern wireless communications is to provide the possibility for high end-user data rates. One straightforward way to meet this requirement is to increase the transmission bandwidth because it is well known that the achievable end-user data rates are limited by transmission bandwidth and transmitter power. However, the transmitted signal in a wireless communication system with wider transmission bandwidth can incur increased corruption due to time dispersion (or equivalently frequency selectivity) of the radio channel. To counteract such signal  1  1.2. Literature Review corruption, researchers have proposed techniques such as receiver-side equalization, multicarrier transmission and specific single-carrier transmission. Orthogonal frequency division multiplexing (OFDM) is a typical multicarrier transmission scheme which can increase the overall transmission bandwidth without suffering from increased signal corruption due to radio channel frequency selectivity. In OFDM systems, data is transmitted in parallel by modulating a number of closely-spaced orthogonal subcarriers, thereby converting a frequency-selective channel into multiple flat fading subchannels. Moreover, intersymbol interference (ISI) can be eliminated by inserting a guard interval between two consecutive OFDM symbols. With these attractive properties OFDM has been adopted by wireless standards such as digital audio broadcasting (DAB), digital video broadcasting (DVB), wireless local area network (WLAN), and wireless metropolitan access network (WMAN) [1–3]. There are several key challenging problems associated with OFDM systems: performance optimization, time and frequency synchronization, channel estimation, and peak-to-average power ratio (PAPR) reduction [4]. In this thesis, we will focus on the channel estimation problem and study channel estimation techniques for an OFDM system operating in dispersive time-varying channels.  1.2 Literature Review Channel estimation is a challenging problem in wireless systems because mobile radio channels are highly dynamic. To avoid channel estimation, one can adopt a differential modulation technique instead of coherent modulation. However, such a system typically results in lower data rates and can incur a 3-4 dB penalty in signal-to-noise ratio (SNR) [5]. For OFDM systems which aim to provide high data rate and spectral efficiency, coherent modulation is more effective; hence, channel estimation is often required as an integral part of the receiver design. Channel estimation for OFDM systems in slow fading channels has been widely studied [6– 8]. However, those channel estimation techniques were developed under the assumption that the channel is constant over one OFDM symbol, an assumption that may not hold in some mobile applications. Channel estimation for OFDM systems in fast fading channels has been studied in 2  1.2. Literature Review the following work. In [9], Li et al. presented a minimum mean square error (MMSE) estimator by exploiting both time-domain and frequency-domain correlations of the frequency response of rapid dispersive fading channels. In a related work, Moon and Choi introduced a channel estimation algorithm by adopting a Gaussian interpolation filter or a cubic spline interpolation filter [10]. However, algorithms in both [9] and [10] require knowledge of channel statistics, which may not be available. To make the estimation algorithm independent of the channel statistics, Li discussed in [11] a robust implementation of the MMSE pilot-symbol-aided estimator, which does not depend on channel statistics. In [12], Zhao and Huang proposed a method employing low-pass filtering in a transform domain to estimate the channel transfer function for the pilot subcarriers. Then a highresolution interpolation is adopted to obtain the channel transfer function for non-pilot subcarriers. In [13], Chang and Su discussed a pilot-aided channel estimation method for OFDM systems operating in Rayleigh fading channels. The channel responses at pilot subcarriers are estimated using a least square (LS) estimator and then interpolated to non-pilot subcarriers using a 2-D regression surface function. This estimator also does not require channel statistics. Recently, Merli and Vitetta proposed a maximum likelihood (ML) based algorithm for joint carrier frequency offset and channel impulse response (CIR) estimation for OFDM systems in [14]. The authors adopted a second-order Taylor series to simplify the estimation problem and derived an approximate closedform solution to the estimation problem. However, it was derived under the assumption that the channel does not change over one OFDM symbol. In [14], the frequency offset is assumed to be constant for all multipaths implying that the Doppler frequency shift is neglected for individual multipath. In this thesis, we study the channel estimation problem for OFDM systems in dispersive timevarying fading channels. Two ML-based channel estimation algorithms which do not require channel statistics are proposed. Unlike the existing solutions in the literature, our channel estimation algorithms are based on a unique channel model which is particularly appropriate for OFDM uplinks in a macrocellular system. In this channel model, the channel state can be determined by estimating a number of unknown channel parameters which do not change over several OFDM symbols. This property of our developed channel model allows us to estimate the channel pa-  3  1.3. Thesis Outline and Contributions rameters once and then employ them to determine the channel state for the next several OFDM symbols.  1.3 Thesis Outline and Contributions This thesis consists of five chapters. Chapter 1 presents background knowledge of wireless communication developments and technologies. In modern wireless communication, a high end-user data rate is one main objective of system design and therefore wider transmission bandwidth is desired. However, an increased transmission bandwidth will increase the frequency selectivity of a radio channel and thus cause severe corruption to the transmitted signal. This problem can be solved by using multicarrier transmission techniques such as OFDM. Chapter 2 provides detailed technical background for the entire thesis. Firstly, fading characteristics of wireless propagation environments are presented and classified into large-scale fading and small-scale fading from the view of propagation terrain. Secondly, the basic concept of OFDM and its implementation are explained. Finally, we list the advantages and drawbacks of OFDM compared with a conventional single-carrier transmission. In Chapter 3, a channel estimation problem is formulated based on a discrete-time CIR model. This channel model is particularly appropriate for OFDM uplinks in a macrocellular system. The channel model developed here allows the CIR to vary within one OFDM symbol. To facilitate the channel estimation problem we adopt a truncated Taylor series expansion to approximate the intercarrier interference (ICI) caused by Doppler frequency shifts. We find that a second-order Taylor series expansion is sufficient for our estimation problem. Then an iterative ML-based algorithm is proposed to estimate the discrete-time channel parameters. Our proposed channel estimation algorithm is particularly useful for an OFDM system which has already compensated the frequency offset due to local oscillator mismatch. The mean square error (MSE) performance of the proposed algorithm in large SNR regions is analyzed using a small perturbation technique, and then demonstrated by simulated results. In Chapter 4, the convergence performance of the proposed algorithm is analyzed. Based on the analysis, an improved fast converging iterative channel estimation algorithm is presented 4  1.3. Thesis Outline and Contributions using a successive overrelaxation (SOR) method to accelerate the proposed algorithm in Chapter 3. Simulated results are also presented to demonstrate the MSE and convergence performance for the improved algorithm. Chapter 5 summarizes the whole thesis and lists our contributions in this work. In addition, some future work related to our current research is suggested.  5  Chapter 2 Wireless Channels and OFDM Technique In this chapter, we will present some background knowledge concerning wireless channels and the OFDM technique. We first address the fading characteristics of wireless propagation environments, then introduce the basic concept and fast Fourier transform (FFT) implementation of OFDM. Finally, advantages and drawbacks of the OFDM technique are compared with conventional single-carrier modulation.  2.1 Wireless Channels In a wireless communication system, the transmitted signal typically undergoes attenuation and distortion over the transmission path. The overall effect on the transmitted signal caused by the transmission path is one main source of system performance degradation in any wireless system. To address and compensate for the attenuation and distortion caused by the transmission path, researchers have studied wireless channels extensively and proposed different channel models [5], [15–18]. These effects, which are mainly due to path loss, shadowing, scattering and reflecting effects caused by unpredictable objects between the transmitter and receiver, can primarily be categorized into large-scale fading and small-scale fading (see Fig. 2.1). In this section, we will briefly describe both types of fading.  2.1.1 Large-Scale Fading Large-scale fading, which is caused by path loss of the signal, can be characterized as a function of transmitted distance and shadowing effects of large obstacles such as buildings and hills. This phenomenon happens when the mobile moves over a large distance (of the order of the cell size),  6  2.1. Wireless Channels  Large-scale fading  Pass loss ) B d( le ve Ll an igS  Small-scale fading  Antenna Displacement Figure 2.1: Signal power fluctuation versus range in wireless channels.  and is typically frequency independent [18]. In an ideal LOS channel with no obstructions between the transmitter and the receiver, the received signal power Pr is given by  Pr = Pt  λ 4π d  2  Gt Gr  (2.1)  where Pt is the transmitted signal power, λ is the wavelength, d is the distance from the transmitter to the receiver, and Gt , Gr are respectively the power gains of the transmit and receive antennas. The power attenuation PL = Pr /Pt is also referred to as free space path loss and it is obvious from (2.1) that the received signal power Pr is inversely proportional to the square of the distance d between the transmit and receive antennas. However, in real transmission environments, the received signal power Pr does not obey this free space path loss model and it varies randomly due to the terrain. Usually a ray tracing method can be used to trace the signal propagation through a wireless channel. Unfortunately, the free space model and ray tracing method cannot model complex propagation environments accurately. Based on empirical measurements, some empirical models for path loss in typical wireless environments have been developed to predict the average received signal power as the transmitted distance d varies, i.e., the Okumura model, the Hata 7  2.1. Wireless Channels model, the European Cooperative for Science and Technology (COST) model, and the piecewise linear (multi-slope) model [17], [19]. In addition to path loss, the transmitted signal is also subject to shadowing, which is caused by changes in reflecting surfaces and scattering objects along the transmission path. The shadowing causes random attenuation to the transmitted signal. Typically, the log-normal shadowing model, which has been confirmed empirically to accurately model the variation in received power, is used to characterize this random attenuation.  2.1.2 Small-Scale Fading Small-scale fading is due to the constructive and destructive addition of different multipath components introduced by the channel between the transmitter and receiver. Therefore, it is also referred to as multipath fading. This phenomenon usually occurs over a distance of several signal wavelengths and is frequency dependent. Since the transmitted signal over a multipath fading channel experiences randomness we must characterize multipath fading channels statistically. Frequency selectivity and the time-varying nature, which depend on the relative relation between parameters of the transmitted signal (i.e., signal bandwidth and symbol duration) and parameters of multipath fading channels (i.e., delay spread and Doppler spread), are two important characteristics of multipath fading channels. In the following, we will briefly describe these two characteristics of multipath fading channels. Depending on the relative relation between transmitted signal bandwidth and delay spread (or equivalently coherence bandwidth BC ), multipath fading channels can be categorized into frequency-nonselective (flat) fading channels and frequency-selective fading channels. The parameter coherence bandwidth is the reciprocal of the delay spread which is defined as the span of the delays of duplicates of the transmitted signal arriving at the receiver via different paths. When the transmitted signal bandwidth is small compared with the coherence bandwidth BC , the channel is called frequency-nonselective or flat fading channel. For a flat fading channel, the spectral components of the transmitted signal are affected in a similar manner so that the multipath components are not resolvable. Otherwise, if the transmitted signal bandwidth is large compared with the coherence bandwidth BC , the channel is said to be frequency-selective. For a frequency-selective fad8  2.2. OFDM Technique ing channel, the spectral components of the transmitted signal are affected by different amplitude gains and phase shifts. In a frequency-selective fading channel, different multipath components with delay differences significantly exceeding the inverse of the transmitted signal bandwidth are resolvable. Typically, such a frequency-selective fading channel can be modeled as a tapped delay line filter with time-variant tap coefficients. In a similar way, the multipath fading channel can be categorized as slow fading or fast fading based on the relative relation between symbol duration and Doppler spread (or equivalently coherence time TC ). The parameter coherence time, which measures the period of time over which the channel effect on the transmitted signal does not change, is defined as the reciprocal of the Doppler spread. The fading channel is said to be slow fading if the symbol duration is small compared with the channel coherence time TC ; otherwise, it is considered to be fast fading. In a slow fading channel, the transmitted signal is affected by the same amplitude gain and phase shift over at least one symbol duration, while the amplitude gain and phase shift vary within one symbol duration in a fast fading channel. According to the discussion above, we have four types of multipath channel, i.e., slow flat fading channel, slow frequency-selective fading channel, fast flat fading channel, and fast frequencyselective fading channel. In this thesis, we focus on the channel state estimation of a fast frequencyselective fading channel for an OFDM system. The channel is characterized using the CIR as [17] L(t)  h(τ ,t) =  ∑ γl (t)δ (τ − τl (t))  (2.2)  l  where L(t) is the number of resolvable multipaths, τl (t) is the delay of the lth multipath, and γl (t) is the corresponding complex amplitude. We will look into the details of this channel model in the next chapter.  2.2 OFDM Technique OFDM has become a promising multicarrier modulation technique and has received considerable interest in the past decade. In OFDM systems, data is transmitted in parallel by modulating a  9  2.2. OFDM Technique number of closely-spaced orthogonal subcarriers, thereby converting a frequency-selective channel into multiple flat fading subchannels [4], [20]. Moreover, ISI can be eliminated by inserting a guard interval between two consecutive OFDM symbols. With these attractive properties OFDM has been adopted by wireless standards such as DAB, DVB, WLAN, and WMAN [1–3]. In this section, we will present the basic concept and FFT implementation of OFDM and list the advantages and drawbacks of the OFDM technique compared with conventional single-carrier modulation.  2.2.1 Basic Concept OFDM is a special case of multicarrier transmission which uses parallel data transmission and frequency division multiplexing. The earliest development of this concept can be traced back to the 1950s [21] and the idea was published in the mid-1960s [22], [23]. In an OFDM system, the available bandwidth is divided into a collection of narrow sub-bands, and the data is transmitted in parallel by modulating these closely-spaced orthogonal subcarriers. Let X[k] denote the complex symbols to be transmitted by an OFDM system, such that the OFDM (modulated) signal can be expressed as N−1  s(t) =  ∑ X[k]e j2π fkt ,  k=0  0 ≤ t ≤ Ts  (2.3)  where fk = f0 + k∆ f , j2 = −1, and N denotes the number of subcarriers in the system. Parameters Ts and ∆ f are called symbol duration and subcarrier spacing of the OFDM system, respectively. These two parameters must satisfy the orthogonality condition (Ts ∆ f = 1) to guarantee that the OFDM signal can be demodulated properly by the receiver. To demonstrate this orthogonal property, we first let  ϕk (t) =  and rewrite (2.3) as     e j2π fk t , 0 ≤ t ≤ Ts    0, otherwise,  k = 0, 1, · · · , N − 1  (2.4)  N−1  s(t) =  ∑ X[k]ϕk (t).  (2.5)  k=0  10  2.2. OFDM Technique Using the definition of ϕk (t), one can show 1 Ts  Ts 0  1 Ts 1 = Ts  ϕk (t)ϕl∗ (t)dt =  Ts 0 Ts  e j2π ( fk − fl )t dt e j2π (k−l)∆ f t dt  0  = δ [k − l]  (2.6)  where δ [k − l] is the delta function defined as  δ [n] =     1, n = 0  (2.7)    0, otherwise.  From (2.6) one can see that {ϕk (t)}N−1 k=0 is a set of orthogonal functions. Using this orthogonal property, the received OFDM signal r(t) (in an ideal case r(t) = s(t)) can be demodulated at the receiver by 1 Ts  Ts 0  r(t)e− j2π fk t dt =  1 Ts  =  1 Ts  Ts 0  s(t)e− j2π fk t dt  Ts 0  N−1  ∑ X[l]ϕl (t)  ϕk∗ (t)dt  l=0  N−1  =  ∑ X[l]δ [l − k]  l=0  = X[k].  (2.8)  According to the discussion above, we provide an illustrative description of a basic OFDM modulator in Fig. 2.2 and the basic principle of OFDM demodulation in Fig. 2.3.  2.2.2 OFDM Implementation Using FFT The relationship between OFDM and the discrete Fourier transform (DFT) was first addressed in [24]. It has been shown that the DFT can be applied to an OFDM system as part of the modulation and demodulation processes. In practice, the DFT can be implemented with the computationally efficient FFT. Here we will briefly discuss the FFT implementation of OFDM. 11  2.2. OFDM Technique  H M πI  W  ;>@  H M πIW  ;>@ ;>@ ;>@  ;>1@  ¦  63  VW  H M πI 1 −W ;>1@  Figure 2.2: Basic principle of OFDM modulation.  H− M πI  W  ³  7V  ³  7V  ³  7V     GW  ;>@  H− M πIW   UW  ;>@  GW  H− M πI 1 −W    GW  ;>1@  Figure 2.3: Basic principle of OFDM demodulation.  12  2.2. OFDM Technique  ;>@  V>@  63  1SRLQW ,))7  ;>1@  ««  ;>1@  V>@  ««  ;>@ ;>@  ;>@  VW 36  '$ &RQYHUWHU  V>1@  Figure 2.4: OFDM modulation by means of IFFT processing.  From (2.3), we can sample s(t) at an interval of Tsa = Ts /N. Then the sampled OFDM signal becomes  s[n] = s(t)|t=nTsa N−1  =  ∑ X[k]e  j2π fk nTs N  ,  k=0  n = 0, 1, . . . , N − 1.  (2.9)  Without loss of generality, we can set f0 = 0 so that fk Ts = k and N−1  s[n] =  ∑ X[k]e  j2π kn N  k=0  = IDFT{X[k]},  n = 0, 1, . . . , N − 1  (2.10)  where IDFT denotes the inverse discrete Fourier transform. Therefore, as illustrated in Fig. 2.4, OFDM modulation can be implemented by means of IDFT (or the computationally efficient inverse fast Fourier transform (IFFT)) processing followed by digital-to-analog conversion. Similarly, the demodulation at the receiver can be carried out using efficient FFT processing, as illustrated in Fig. 2.5.  2.2.3 Advantages and Drawbacks of OFDM As a special case of multicarrier transmission scheme, OFDM has the following key advantages over conventional single-carrier transmission [20]:  13  2.2. OFDM Technique  Q7VD  63  U>@  ;>@  U>1@  1SRLQW ))7  ««  U>Q@  ;>@  ««  UW  U>@  ;>@ ;>@  ;>1@  36  ;>1@  Figure 2.5: OFDM demodulation by means of FFT processing.  • OFDM is an efficient way to deal with multipath fading; for a given channel delay spread, the implementation complexity is much lower than that of a conventional single-carrier system with an equalizer. • In relatively slow time-varying channels, it is possible to significantly enhance the capacity by adapting the data rate per subcarrier according to the SNR of that particular subcarrier. • OFDM is robust against narrowband interference, because such interference only affects a small percentage of the subcarriers. • OFDM makes single-frequency networks possible, which is especially attractive for broadcasting applications. However, OFDM also has some drawbacks compared with conventional single-carrier modulation [20]: • OFDM is more sensitive to frequency offset and phase noise. The frequency offset and phase noise will destroy the orthogonality among subcarriers and hence introduces ICI. • OFDM signals have relatively large PAPR, which tends to reduce the power efficiency of the RF amplifier. • OFDM needs an adaptive or coded scheme to overcome spectral nulls in the channel.  14  Chapter 3 ML-based Channel Estimation for OFDM Systems in Dispersive Time-Varying Channels In this chapter, by taking into account the ICI component, we propose an ML-based channel estimation algorithm (which does not require channel statistics) for OFDM systems in dispersive time-varying fading channels. A channel estimation problem is formulated based on a discretetime CIR model. This channel model is particularly appropriate for OFDM uplinks in a macrocellular system. The channel model developed here allows the CIR to vary within one OFDM symbol. To facilitate the channel estimation problem we adopt a truncated Taylor series expansion to approximate the ICI caused by Doppler frequency shifts. We find that a second-order Taylor series expansion is sufficient for our estimation problem. Then an iterative ML-based algorithm is proposed to estimate the discrete-time channel parameters. Using a small perturbation technique, we analyze the performance of the proposed algorithm in large SNR regions. Our computer simulations demonstrate that the proposed algorithm can estimate the desired parameters accurately, and the simulated performances agree with our theoretical analysis.  3.1 System Model We consider an OFDM system with N subcarriers. For each OFDM symbol, we denote the transmitted symbols as X[0], X[1], . . . , X[N − 1]. After the IFFT, the time-domain OFDM signal can be  15  3.2. Formulation of Channel Estimation Problem expressed as [25] s[n] =  j2π kn 1 N−1 X[k]e N , ∑ N k=0  n = 0, 1, . . . , N − 1  (3.1)  where j2 = −1. Each OFDM symbol is extended with a cyclic prefix and transmitted after appropriate pulse shaping. Recalling (2.2), we can describe a multipath CIR by L  h (τ ,t) =  ∑ γl (t) δ (τ − τl )  (3.2)  l=1  where γl (t) is the complex amplitude of the lth multipath. The number of resolvable multipaths L(t) and the delay of the lth multipath τl (t) are assumed to be L and τl , respectively. According to [6], [9], the CIR can be well approximated by a sampled discrete-time CIR h[m, n] = h(mTsa , nTsa ), where Tsa is the sampling interval defined as Tsa = Ts /N and Ts is the OFDM symbol duration. In this thesis, we assume a time-varying multipath channel and express the discrete-time CIR as [12], [26] L  h[m, n] =  ∑ hl e l=1  j2π fl n N  δ [m − nl ]  (3.3)  where hl is the complex amplitude for the lth multipath, fl is the corresponding Doppler frequency shift normalized by 1/Ts , and nl is the delay in samples for the lth multipath. The detailed derivation of this discrete-time CIR model is given in Appendix A. It should be emphasized that the discrete-time CIR in (3.3) is only valid for OFDM uplink transmission in macrocellular systems. Without loss of generality, we assume n1 ≤ n2 ≤ · · · ≤ nL . It is seen from (3.3) that the CIR can vary within one OFDM symbol duration.  3.2 Formulation of Channel Estimation Problem Our goal is to estimate the channel parameters hl , fl , for l = 1, 2, . . . , L, in the discrete-time CIR in (3.3). We assume that the length of the cyclic prefix is longer than the delay spread of the multipath channel and that the receiver accurately removes the cyclic prefix before implementing the FFT.  16  3.2. Formulation of Channel Estimation Problem For the channel model given in (3.3), one can express the output of the FFT as  Y [k] =  j2π k (n−nl ) j2π fl n j2π kn 1 L N−1 N−1 X[k ]e N hl e N e− N + w[k] ∑ ∑ ∑ N l=1 k =0 n=0  L N−1  =  ∑∑  e−  j2π k nl N  l=1 k =0 L N−1  =  1 − e j2π (k −k+ fl )  X[k ]hl  N 1−e  ∑ ∑ e−  j2π k nl N  j2π (k −k+ fl ) N  + w[k]  X[k ]hl R(k , k, fl ) + w[k]  l=1 k =0  = αk X[k] + βk + w[k],  k = 0, 1, . . . , N − 1  (3.4a)  where L  αk =  ∑ e−  j2π knl N  hl R (k, k, fl )  (3.4b)  l=1  and L N−1  βk =  ∑ ∑ e−  j2π k nl N  X[k ]hl R(k , k, fl ).  l=1 k =k k =0  (3.4c)  In (3.4a) w[k] is an additive complex Gaussian noise component with zero mean and variance σ 2 , i.e., CN 0, σ 2 , and  R(k , k, fl ) = =  1 N−1 j2π n(k −k+ fl ) ∑e N N n=0 1 − e j2π (k −k+ fl ) N 1−e  j2π (k −k+ fl ) N  (3.5)  represents the interference of X[k ] on the kth subcarrier caused by a normalized Doppler frequency shift fl . In (3.4b) and (3.4c), αk and βk denote the multiplicative distortion at the desired subchannel and the additive ICI term, respectively. In the absence of Doppler frequency shifts ( fl = 0, l = 1, 2, . . . , L), we have R(k , k, 0) = δ [k − k] resulting in the conventional model without ICI among subcarriers. In the special case when fl = ε , for l = 1, 2, . . . , L, our estimation problem  17  3.2. Formulation of Channel Estimation Problem is mathematically equivalent to the joint channel and frequency offset estimation problem that has been considered in [14], [27] and [28]. In practical mobile communication scenarios, fl is typically less than 0.1. Consider an OFDM system, for example, with N = 256 subcarriers and a subcarrier spacing of 7.81 kHz. A mobile terminal moving at a speed of 84.4 km/h will result in a normalized Doppler frequency fl = 0.05 with a carrier frequency fc = 5 GHz. Exploiting the fact that fl has small values, one can consider a Taylor series expansion for R(k , k, fl ) as  R(k , k, fl ) =   2   1 + jπ N−1 fl − π 2 2N 2 −3N+1 fl2 − jπ 3 (N−1) fl3 + · · · , N 3N 2 3N 2  k =k  2 2 +3ω (1+ω −N ω +N)   − j2π fl + 2π 2 (2−N)ω2+N f 2 + j4π 3 N (1−ω ) 3N fl3 + · · · , k = k 3 (1−ω )3 (1−ω )N (1−ω ) N 2 l (3.6a)  where  ω = exp  j2π (k − k) . N  (3.6b)  Applying (3.6a) to (3.4a) and collecting the terms with the same power of fl , one obtains L  Y [k] =  Al [k] + Bl [k] fl +Cl [k] fl2 + Dl [k] fl3 + · · · + w[k]  ∑ hl l=1  (3.7)  for k = 0, 1, . . . , N − 1, or, in vector form, L  Y=  ∑ hl l=1  Al + Bl fl + Cl fl2 + Dl fl3 + · · · + w.  (3.8)  The first three leading coefficients Al [k], Bl [k], Cl [k] in (3.7) can be shown to be  Al [k] = X[k]e  − j2π knl N  (3.9a)  and Bl [k] = jπ  − j2π k nl − j2π knl N −1 j2π N−1 X[k ]e N X[k]e N − ∑ N (1 − ω )N k =k  (3.9b)  k =0  18  3.2. Formulation of Channel Estimation Problem and  Cl [k] = −π 2  N−1 − j2π k nl − j2π knl 2N 2 − 3N + 1 2 (2 − N)ω + N N X[k ]e N . X[k]e π + 2 ∑ 2 2 2 3N (1 − ω ) N k =k  (3.9c)  k =0  The channel estimation technique in this thesis will be based on (3.8). In (3.8), the coefficient vectors Al , Bl , Cl , . . . can be derived given the transmitted training sequence X[k], k = 0, 1, . . . , N − 1. We will assume the receiver has knowledge of the training sequence X[k] and performs parameter estimation of the hl ’s and fl ’s based on (3.8). We note that the Taylor expansion in (3.6a) is an approximation for R(k , k, fl ) and the accuracy of the approximation depends on the order of the Taylor expansion. On the other hand, the order of the Taylor expansion also determines the computational complexity of the channel estimation algorithm. In this work we employ a second-order Taylor expansion for | fl | < 0.1. In order to assess the accuracy of the second-order Taylor series expansion, we compare the MSE of the received OFDM symbol using several approximation orders. Here the MSE is defined as the mean squared difference between the exact value of the signal at the receiver (Y [k]) and its approximation (Y [k])  MSE = E  Y [k] − Y [k]  2  (3.10)  where E{x} is the ensemble average of a random sequence x. Y [k] and Y [k] are calculated according to (3.4a) and (3.7) assuming w[k] = 0, respectively. From Fig. 3.1 one observes that the MSE increases with the normalized Doppler frequency fl as expected. Fig. 3.1 also shows that the first-order approximation may be inadequate. However, the difference between the second-order approximation and third-order approximation is negligible. Furthermore, with a typical value of fl , e.g. fl = 0.05, the MSE of the second-order approximation is about −50 dB and, thus, can be ignored in our estimation problem. Therefore, we can conclude that a second-order approximation  19  3.2. Formulation of Channel Estimation Problem  0 −10 −20  MSE (dB)  −30 −40 −50 −60 −70 −80 p=1 p=2 p=3  −90 −100  0.02  0.04  0.06 0.08 0.1 0.12 fl (Normalized Doppler Frequency)  0.14  Figure 3.1: An MSE comparison for three different approximation orders.  is sufficient for our problem and modify (3.8) to be L  Y≈  ∑ hl  Al + Bl fl + Cl fl2 + w.  (3.11)  l=1  The likelihood function of the received symbol is then given by  f Y  {hl , fl }Ll=1    L 1 = N N exp − Y − ∑ hl Al + Bl fl + Cl fl2  π σ l=1 L  H  × Y − ∑ hl Al + Bl fl + Cl fl2  (3.12) /σ 2  l=1  20  3.2. Formulation of Channel Estimation Problem where (·)H denotes the Hermitian transpose of the argument matrix. One can therefore express the ML estimation for the parameters hl ’s and fl ’s as L hˆ l , fˆl l=1  2  L  = arg min  {hl , fl }Ll=1  Y − ∑ hl  Al + Bl fl + Cl fl2  (3.13)  l=1  where hˆ l ’s and fˆl ’s are the solutions to  ∂ ln f Y {hl , fl }Ll=1 = 0 ∂ hl  (3.14)  ∂ ln f Y {hl , fl }Ll=1 = 0, ∂ fl  (3.15)  and  respectively. From (3.12), we can get  ln f Y  L 1 1 = ln N N − 2 Y − ∑ hl Al + Bl fl + Cl fl2 π σ σ l=1  {hl , fl }Ll=1  H  (3.16)  L  × Y − ∑ hl Al + Bl fl + Cl fl2 l=1  and show that (3.14) and (3.15) are equivalent to the following equations H  L  Y − ∑ hl  Al + Bl fl + Cl fl2  Al + Bl fl + Cl fl2 = 0  (3.17)  l=1  ℜ      L  Y − ∑ hl Al + Bl fl + Cl fl2 l=1  H  hl Bl + 2Cl fl      = 0,  (3.18)  respectively. Here ℜ{·} denotes the real part of the complex argument. It is worth mentioning that we need to pay attention to the properties of complex derivatives [29] when obtaining (3.17) and (3.18). Furthermore, we also need to note that (3.17) and (3.18) are non-linear in the hl ’s and fl ’s. A technique to solve these two equations will be presented in the next section.  21  3.3. An Iterative Channel Estimation Algorithm Table 3.1: Algorithm 1 (0)  Initialization: Initialize hl  (0)  and fl , for l = 1, 2, . . . , L. Set the iteration counter i = 0.  Iterations: Step 1: For each l = 1, 2, . . . , L, calculate (i) (i+1) h˜ l = arg min L hl , f˜l hl  (i+1) f˜l  (i+1) = arg min L h˜ l , fl  Step 2: If max  fl ˜h(i+1) −h˜ (i) l l (i) h˜ l  l  (i+1) ˜(i) f˜l − fl  × 100% > δh or max  (i) f˜l  l  (i+1)  i = i + 1 and go to Step 1, otherwise output h˜ l  × 100% > δ f , let  (i+1) ’s and f˜l ’s.  3.3 An Iterative Channel Estimation Algorithm In this section, we propose an iterative approach to estimate the hˆ l ’s and fˆl ’s based on (3.13). The most straightforward approach is the coordinate-ascent algorithm [30] as shown in Table 3.1, where the cost function L(hl , fl ) is defined as l−1  (i+1)  L (hl , fl ) = Y − ∑ h j  (i+1)  Aj +Bj fj  (i+1) 2  +Cj f j  j=1 L  −  ∑ j=l+1  (i)  hj  (i)  Aj + Bj f j +Cj  (i) 2 fj  (3.19)  2  − hl Al + Bl fl +Cl fl 2  .  Note that in Algorithm 1, the first equation in Step 1 is a linear least-square problem and the second equation in Step 1 requires finding the roots for a third-order polynomial. Both equations can be solved using well-known numerical methods. As discussed in [30], the coordinate-ascent algorithm is a special case of the space-alternating generalized expectation-maximization algorithm and has the attractive property of increasing the likelihood value at each iteration.  22  3.4. Performance Analysis of the Channel Estimation Algorithm  3.4 Performance Analysis of the Channel Estimation Algorithm In this section, we analyze the MSE performance of the proposed algorithm in large SNR regions assuming the fl ’s are small. From the discussion in Section 3.2 we can see that the ML-based estimation algorithm aims to minimize the cost function in (3.13). Using a small perturbation technique [31], one can obtain from (3.13) a linear equation for the perturbation of the estimates as    R11 R12   R21 R22   . ..  .. .   RL1 RL2 where    R  ℑ   ℜ {∆hl }  ∆θ l =   ℑ {∆hl }  ∆ fl  AHj Ai  ℜ AHj Ai ℑ hi AHj Bi          (3.20a)  n  ∆θ  AHj Ai        · · · R1L   ∆θ1   n1          · · · R2L    ∆θ2   n2   . = .  ..  ..     . .    ..   ..      nL · · · RLL ∆θ L   ℜ  H Ri j =   −ℑ A j Ai  ℜ hi AHj Bi  and      ℜ  h∗j BHj Ai  −ℑ h∗j BHj Ai ℜ{hi h∗j BHj Bi } wH Al   ℜ  H nl =   −ℑ w Al  ℜ hl wH Bl                (3.20b)  (3.20c)  where ℜ{·} and ℑ{·} respectively denote the real part and imaginary part of the argument, {·}∗ denotes the complex conjugate of the argument, ∆hl = hˆ l − hl and ∆ fl 1 = fˆl − fl are the estimate perturbations caused by Gaussian noise nl . A detailed derivation of (3.20) is given in Appendix B. When the Gaussian noise is absent, we have hˆ l = hl and fˆl = fl . Equation (3.20) indicates that the perturbations ∆θl subject to the Gaussian noise components 1 Here  ∆ fl does not denote the subcarrier spacing.  23  3.4. Performance Analysis of the Channel Estimation Algorithm nl , l = 1, 2, . . . , L, are also Gaussian, and one can express ∆θl as        ∆θ1   R11 R12     ∆θ2   R21 R22     . = . ..  ..   .. .       ∆θ L RL1 RL2  −1   · · · R1L   · · · R2L   ..  .. . .    · · · RLL             n1   n2   ..  .    nL  (3.21)  or in a more compact form ∆θ = R−1 n.  (3.22)  To determine the variance of ∆θl subject to Gaussian noise components nl , we have  H HA HA ℜ w ℜ w i j       H    H H ni n j =  −ℑ w Ai   −ℑ w A j      ℜ hi wH Bi ℜ h j wH B j  H H −ℜ wH Ai ℑ wH A j  ℜ w Ai ℜ w A j  H H = ℑ wH Ai ℑ wH A j  −ℑ w Ai ℜ w A j  −ℜ hi wH Bi ℑ wH A j ℜ hi wH Bi ℜ wH A j     ℜ wH Ai ℜ h j wH B j −ℑ  wH Ai  ℜ hj  wH B  j  ℜ hi wH Bi ℜ h j wH B j        (3.23)  and  E ni nHj    H     ℜ A j Ai  2  σ H E  =  −ℑ A j Ai 2       ℜ hi AH Bi j  σ2 E Ri j 2 σ2 = Ri j . 2 =  ℑ  AHj Ai  ℜ  h∗j BHj Ai  ℜ AHj Ai  −ℑ h∗j BHj Ai  ℑ hi AHj Bi  ℜ hi h∗j BHj Bi              (3.24)  In obtaining (3.24), we have used the fact that w is a complex Gaussian noise vector whose elements have independent real and imaginary parts with zero mean and variance equal to σ 2 /2.  24  3.4. Performance Analysis of the Channel Estimation Algorithm Therefore, using (3.20a), (3.22) and (3.24), one can show straightforwardly that        E          and     ∆θ 1     . . E   .       ∆θL  n1 .. .             H σ2  n1 · · · nH R = L   2      nL  (3.25)             σ 2 −1  ∆θ1H · · · ∆θLH R . =   2       (3.26)  Generally speaking, the elements in Ri j , i = j, are non-zero and, thus, the small perturbations for different taps, ∆θl , are correlated. When L << N and a pseudo-random sequence X[k] is used as the training symbols, the values of Ri j , i = j, are generally much smaller than those in Rii and, thus, negligible [25]. For such cases, one can rewrite (3.26) as  E ∆θlH ∆θl =  σ 2 −1 R , 2 ll  l = 1, 2, . . . , L  (3.27)  where   AH l Al   ℜ  H Rll =   −ℑ Al Al  ℜ hl AH l Bl  AH l Al   = 0   ℜ hl AH l Bl  ℑ  AH l Al  ℜ AH l Al ℑ hl AH l Bl 0 AH l Al ℑ hl AH l Bl  ℜ  h∗l BH l Al  −ℑ h∗l BH l Al ℜ{h∗l hl BH l Bl } ℜ h∗l BH l Al −ℑ  h∗l BH l Al  |hl |2 BH l Bl }                (3.28)  25  3.4. Performance Analysis of the Channel Estimation Algorithm and R−1 ll is shown to be R−1 ll =  1 2  H |hl |2 AH l Al ζ − Bl Al  2 2 HB | |h ζ − ℑ h A l l l l   H H ×  ℜ hl Al Bl ℑ hl Al Bl  H −AH l Al ℜ hl Al Bl  H ℜ hl AH l Bl ℑ hl Al Bl  ζ |hl |2 − ℜ hl AH l Bl  2  H −AH l Al ℜ hl Al Bl H −AH l Al ℑ hl Al Bl  H −AH l Al ℑ hl Al Bl  AH l Al  2         (3.29)  H where ζ is defined as ζ = AH l Al Bl Bl for the sake of convenience. Equation (3.27) indicates that  the perturbations ∆θl for different taps can be treated independently. However, for the same l, the perturbations ∆hl = ℜ{∆hl } + jℑ{∆hl } and ∆ fl are still correlated because the non-diagonal elements in Rll are not negligible. It is of interest to obtain some explicit results for the estimation performance. Motivated by [32], we consider the asymptotic behavior with a white training sequence. Under this condition, we have AH l Al = 1 BH l Al = − BH l Bl =  (3.30a)  jπ (N − 1) N  (3.30b)  2π 2 (N − 1)(2N − 1) 3N 2  (3.30c)  where the derivations of (3.30a)-(3.30b) are trivial, and a detailed derivation of (3.30c) can be found in Appendix C. According to the definition of ∆θl , we can obtain     H   ℜ {∆hl }   ℜ {∆hl }       ∆θl ∆θlH =   ℑ {∆hl }   ℑ {∆hl }     ∆ fl ∆ fl  2 ℜ {∆hl } ℑ {∆hl } ℜ {∆hl } ∆ fl  (ℜ {∆hl })  = ℑ {∆hl } ∆ fl (ℑ {∆hl })2  ℜ {∆hl } ℑ {∆hl }  ∆ fl2 ℑ {∆hl } ∆ fl ℜ {∆hl } ∆ fl         (3.31)  26  3.4. Performance Analysis of the Channel Estimation Algorithm and ∆h∗l ∆hl = |∆hl |2 = (ℜ {∆hl })2 + (ℑ {∆hl })2 .  (3.32)  Therefore, we can obtain the expression for the estimation performance as E {∆h∗l ∆hl } = E (ℜ {∆hl })2 + (ℑ {∆hl })2 = E (ℜ {∆hl })2 + E (ℑ {∆hl })2 .  (3.33)  Applying (3.27), (3.29) and (3.30) to (3.33), one can obtain  E {∆h∗l ∆hl } = =  H H σ 2 2AH l Al Bl Bl − Bl Al  2  H H H 2AH l Al Al Al Bl Bl − Bl Al  2  (5N 2 − 6N − 1)σ 2 . 2(N 2 − 1)  (3.34)  Similarly, we can obtain  E ∆ fl2 =  =  σ 2 AH l Al H H 2 |hl |2 AH l Al Bl Bl − Bl Al  3N 2 σ 2 2π 2 (N 2 − 1) |hl |2  2  .  (3.35)  When N → ∞, we have (5N 2 − 6N − 1)σ 2 5σ 2 = N→∞ 2(N 2 − 1) 2  lim E {∆h∗l ∆hl } = lim  N→∞  (3.36)  and  lim E ∆ fl2 = lim  N→∞  3N 2 σ 2  N→∞ 2π 2 (N 2 − 1) |hl |2  =  3σ 2 2π 2 |hl |2  .  (3.37)  27  3.5. Numerical Results and Discussions In comparison, the Cramer-Rao bound (CRB) at large SNR for the classical frequency offset estimator of Moose is [33], [34]  CRBMoose ( fl ) =  2 2π  2  2σ 2 |hl |2 AH l Al  =  2σ 2  π 2 |hl |2  .  (3.38)  Thus, from (3.37) and (3.38), we can see that the estimation algorithm considered in this paper can 2 = 1.25 dB gain over Moose’s estimator. In Section 3.5, it will be seen that achieve 10 log10 3/2  simulations performed to test the analytical results confirm the theoretical analysis.  3.5 Numerical Results and Discussions In this section, we present some numerical results by considering a discrete baseband OFDM uplink in our simulation. The complex amplitude hl is randomly chosen, and the normalized Doppler frequency fl is chosen to be less than 0.1. The total number of subcarriers is fixed at N = 256. A binary phase-shift keying (BPSK) signaling scheme is adopted in the simulation. We set the stopping threshold values for the proposed algorithm as δh = 1% and δ f = 1%. A second-order Taylor series expansion is adopted in the simulation. Here the MSE is adopted as the performance measure. Fig. 3.2 plots the MSE performance comparison between Algorithm 1 and Moose’s estimator for fl in a fading channel with L = 1. Here hl is randomly chosen and fl is set to be fl = 0.05. The MSEs are obtained by averaging over 2000 trials, and the theoretical performance predictions of fˆl are obtained using (3.37) and (3.38). From Fig. 3.2, we observe that the simulated MSE performance of both algorithms and their theoretical predictions are in excellent agreement over a wide range of SNR values. The simulated performance curve for Algorithm 1 exhibits an error floor, which is caused by the approximation error when the average SNR value is greater than 25 dB. This is the cause of discrepancy between the simulated curve and the theoretical curve for SNR greater than 25 dB. However, the values of MSE in these regions are small and the error floor is negligible. Fig. 3.2 also shows that Algorithm 1 outperforms Moose’s estimator by about 1.2 dB over a wide range of SNR values. 28  3.5. Numerical Results and Discussions Fig. 3.3 plots the average MSE performance of hˆ l in a multipath fading channel with L = 3 for the proposed algorithm and Figs. 3.4-3.6 plot MSE performance of fˆl for the same system with three different fl values. Here the hl ’s are arbitrarily chosen to be hl = [0.2944 + 1.6236 j, −1.3362 − 0.6918 j, 0.7143 + 0.858 j], and the fl ’s are set to be fl = [0.04, 0.02, 0.06], respectively. The MSEs are obtained by averaging over 1500 trials, and the theoretical performance predictions of hˆ l and fˆl are obtained using (3.36) and (3.37), respectively. From Figs. 3.3-3.6, we observe that the simulated MSE performance of hˆ l and fˆl for Algorithm 1 and their theoretical predictions have excellent agreement over a wide range of SNR values. Figs. 3.4-3.6 also show that the MSE performance degrades with an increasing value of fl , owing to increased (deterministic) approximation errors. From Fig. 3.4, it is seen for fl = 0.02 that the simulated result and the theoretical result agree over the entire SNR range. Figs. 3.3, 3.5 and 3.6 suggest that the simulated results overestimate the theoretical results when the average SNR value is large. This is because we adopt a second-order Taylor series expansion to simplify the estimation problem, and the error floor is caused by the approximation error. However, the values of MSE in these regions are small and the error floor is negligible.  29  3.5. Numerical Results and Discussions  −10 Simulation result for Algorithm 1 Simulation result for Moose estimator Theoretical value for Algorithm 1 Theoretical value for Moose estimator  −20  MSE of fˆl (dB)  −30  −40  −50  −60  −70  −80 −10  0  10 20 Average SNR (dB)  30  40  Figure 3.2: An MSE performance comparison between Algorithm 1 and Moose’s estimator for fˆl when fl = 0.05, δh = 1% and δ f = 1%.  30  3.5. Numerical Results and Discussions  0 Simulated result Theoretical value  ˆl (dB) Average MSE of h  −10  −20  −30  −40  −50  −60 −10  0  10 20 Average SNR (dB)  30  40  Figure 3.3: An average MSE performance comparison between the simulated results and the theoretical values for hˆ l of Algorithm 1 when δh = 1% and δ f = 1%.  31  3.5. Numerical Results and Discussions  −10 Simulated result Theoretical value −20  MSE of fˆl (dB)  −30  −40  −50  −60  −70 −10  0  10 20 Average SNR (dB)  30  40  Figure 3.4: An MSE performance comparison between the simulated results and the theoretical values for fˆl of Algorithm 1 when fl = 0.02, δh = 1% and δ f = 1%.  32  3.5. Numerical Results and Discussions  −10 Simulated result Theoretical value −20  MSE of fˆl (dB)  −30  −40  −50  −60  −70 −10  0  10 20 Average SNR (dB)  30  40  Figure 3.5: An MSE performance comparison between the simulated results and the theoretical values for fˆl of Algorithm 1 when fl = 0.04, δh = 1% and δ f = 1%.  33  3.5. Numerical Results and Discussions  −10 Simulated result Theoretical value −20  MSE of fˆl (dB)  −30  −40  −50  −60  −70 −10  0  10 20 Average SNR (dB)  30  40  Figure 3.6: An MSE performance comparison between the simulated results and the theoretical values for fˆl of Algorithm 1 when fl = 0.06, δh = 1% and δ f = 1%.  34  Chapter 4 Convergence Study of ML-based Channel Estimation Algorithms In Chapter 3, by taking into account the ICI component, we proposed an ML-based channel estimation algorithm (which does not require channel statistics) for OFDM systems in dispersive timevarying fading channels. Our proposed channel estimation algorithm is iterative. In this chapter, we will study the convergence performance of the iterative algorithm presented in Chapter 3 analytically. Furthermore, we propose an improved fast converging channel estimation algorithm for OFDM systems in dispersive time-varying fading channels. An SOR method [35] is adopted to accelerate our proposed channel estimation algorithm. Our computer simulations demonstrate that the improved algorithm, which achieves the same MSE performance as the algorithm proposed in Chapter 3, has better convergence performance.  4.1 Convergence Performance Analysis of the Proposed Algorithm In this section we first analyze the convergence rate of Algorithm 1 proposed in Chapter 3. Based on the convergence rate analysis we propose an improved fast converging algorithm using an SOR method [35] in the next section. We first rewrite the linear expression in (3.20a) for small pertur-  35  4.1. Convergence Performance Analysis of the Proposed Algorithm bation as a linear Gauss-Seidel iteration [35, eq. (3.15), page 72]     (i+1)   ∆θ 1   ∆θ (i+1)  2 D ..  .   (i+1) ∆θL    (i+1)    ∆θ1     ∆θ (i+1)   2  = −L  ..   .     (i+1) ∆θL      (i) θ ∆ 1           ∆θ (i)     2    −U  .  +    ..           (i) ∆θL    n1   n2   ..  .    nL  (4.1)  or, equivalently, D∆θ (i+1) = −L∆θ (i+1) −U∆θ (i) + n  (4.2)  where D is a diagonal matrix whose entries are the diagonal elements of R; L and U are the lower and upper triangular matrices of R, i.e., L + D +U = R. Applying (3.20a) to (4.1), we can obtain   (i+1)   ∆θ 1   ∆θ (i+1)  2 D ..  .   (i+1) ∆θ L i.e.,      (i+1)    ∆θ1     ∆θ (i+1)   2  = −L  ..   .     (i+1) ∆θL   (i+1)  − ∆θ1  ∆θ 1   ∆θ (i+1) − ∆θ2  2 (D + L)  ..  .   (i+1) − ∆θL ∆θL  or, in a more compact form as       (i) θ θ ∆ ∆ 1  1          (i)   ∆θ   ∆θ2      2   −U  .  + R  .  ,   .   .    .   .       (i) ∆θ L ∆θL       (i)   ∆θ1 − ∆θ1     ∆θ (i) − ∆θ2    2  = −U  ..   .     (i) ∆θL − ∆θL  (4.3)            (D + L) ∆θ (i+1) − ∆θ = −U ∆θ (i) − ∆θ .  (4.4)  (4.5)  The convergence rate2 of (4.5) is determined by the largest (in an absolute value sense) eigenvalue of (D + L)−1U [35]. The Gauss-Seidel iteration is a simple and efficient algorithm. However, by exploiting the 2 With  this convention, the smaller the largest eigenvalue is, the faster the iterative algorithm will converge.  36  4.1. Convergence Performance Analysis of the Proposed Algorithm structure of the Fisher information matrix R, one may adopt an overrelaxation philosophy to accelerate the convergence rate. We will next introduce an SOR method [35] to achieve faster convergence. Corresponding to the linear iteration (4.2), an overrelaxed linear iteration can be expressed as [35, eq. (3.22), page 73] D∆θ (i+1) = η −L∆θ (i+1) −U∆θ (i) + n + (1 − η ) D∆θ (i)  (4.6)  where η is the relaxation factor. Applying (3.20a) to (4.6), we can obtain (D + η L) ∆θ (i+1) − ∆θ = [(1 − η )D − η U] ∆θ (i) − ∆θ .  (4.7)  Similarly to the Gauss-Seidel iteration in (4.5), the convergence rate of the successive overrelaxed iteration in (4.7) depends on the largest absolute value of the eigenvalues of (D + η L)−1 [(1 −  η )D − η U] [35]. To achieve the fastest convergence rate for (4.7), one is required to minimize the largest eigenvalue of (D + η L)−1 [(1 − η )D − η U] by choosing the optimal value for η . To simplify the analysis, we again assume that Ri j = 0 for i = j. Under this assumption, D, L and U are block diagonal matrices composed of 3 × 3 matrices Dl , Ll and Ul , l = 1, . . . , L. Thus, for the first step, we consider the optimization on the eigenvalues of (Dl + η Ll )−1 [(1 − η )Dl − η Ul ] separately. For the sake of convenience, we denote al = AH l Al , bl = |hl |2 BH l Bl , cl = ℜ h∗l BH l Al , dl = ℑ h∗l BH l Al ,  37  4.1. Convergence Performance Analysis of the Proposed Algorithm and obtain  Mη = Dl + η Ll   0 0   al    = 0 a 0 l     η cl −η dl bl   1 0 0   al 0 0     0 a 0 = 0 1 0 l     0 0 bl η cl /al −η dl /al 1  Nη = (1 − η )Dl − η Ul  0 −η cl  (1 − η )al  = η dl 0 (1 − η )al   0 0 (1 − η )bl         (4.8)      .    (4.9)  It is straightforward to show that  Mη−1 Nη      =    1−η  0  −η cl /al  0  1−η  η dl /al  −η (1 − η )cl /bl η (1 − η )dl /bl (1 − η ) + η 2 (c2l + dl2 )/al bl      .    (4.10)  The eigenvalues of (4.10) are the roots of the characteristic equation  λ 3 − 3(1 − η ) +  η 2 (c2l + dl2 ) 2 η 2 (c2l + dl2 ) λ + (1 − η ) 3(1 − η ) + λ − (1 − η )3 = 0. (4.11) al bl al bl  For the sake of simplicity, we denote 2  BH l Al c2l + dl2 = H kl = al bl Al Al BH l Bl  (4.12)  38  4.1. Convergence Performance Analysis of the Proposed Algorithm and simplify (4.11) to  λ 3 − 3(1 − η ) + η 2 kl λ 2 + (1 − η ) 3(1 − η ) + η 2 kl λ − (1 − η )3 = 0.  (4.13)  It is shown in Appendix D that the optimal value of η which minimizes the largest eigenvalue (in absolute value) of Mη−1 Nη , is  √ 2(1 − 1 − kl ) η= , kl  (4.14)  and the corresponding largest absolute value for the eigenvalues of Mη−1 Nη is √ 2(1 − 1 − kl ) . |λ | = |1 − η | = 1 − kl  (4.15)  Generally, an SOR method requires knowledge about the structure of the Fisher information matrix. It is observed in (3.20b) that this matrix (the Fisher information matrix Rll ) depends on the channel realization hl . However, (4.12) and (4.14) indicate that for a given 3 × 3 Fisher information matrix Rll , the optimal η value does not depend on the channel realization hl . Instead, it only relies on the training sequence. In other words, given the training sequence, the receiver can employ η in (4.14) to achieve a fast convergence rate for the overrelaxed iteration. Typically convergence accelerating algorithms require full knowledge of the Fisher information matrix Rll to achieve better convergence rates. This implies that one needs to estimate hl and it renders such algorithms less useful. In contrast, the overrelaxation method proposed in this work only requires information about the training sequence and, therefore, is more robust. This property makes the overrelaxation method an attractive approach to achieve a fast convergence rate in practical receivers. From (4.12) and (4.14), it is observed that the η value depends on kl and may vary for different l. Substitution of (3.30) into (4.12) yields  π 2 (N − 1)2 /N 2 2π 2 (N − 1)(2N − 1)/3N 2 3(N − 1) = 2(2N − 1)  kl =  (4.16)  39  4.1. Convergence Performance Analysis of the Proposed Algorithm which is, in fact, independent of l assuming a white training sequence is employed. By letting N → ∞, we have the limiting value η , from (4.14) and (4.16), as √ 2 1 − 1 − kl η = lim N→∞ kl 4 = 3  (4.17)  which leads to a convergence rate of |1 − η | = 1/3. In comparison, considering the Gauss-Seidel iteration by letting η = 1, one can express (4.10) as     0 0 0     −1  Mη Nη =  0 0 0     0 0 kl  (4.18)  and conclude that the largest eigenvalue for the iteration is kl . Thus, the largest eigenvalue for the Gauss-Seidel iteration is kl = 3/4, which is larger than that of 1/3 for the SOR method, implying a faster convergence for the later3 . In each Gauss-Seidel iteration, a value that maximizes the likelihood function is produced as an output. In contrast, an SOR method adopts a less aggressive strategy by combining the GaussSeidel solution and the previous value linearly as (i+1)  (i)  (i+1)  = η h˜ l  (i+1)  (i) (i+1) + (1 − η ) fl = η f˜l  hl  fl  + (1 − η )hl  (4.19a)  (4.19b)  (i+1) (i+1) and f˜l are the solutions to the Gauss-Seidel where η is an overrelaxation parameter, and h˜ l  iteration. When η = 1, the SOR method degenerates to a Gauss-Seidel iteration. 3 We  should comment that the η value in (4.17) is only asymptotically optimal. In Section 4.3, we will demonstrate that this suboptimal solution already shows significant improvement over the Gauss-Seidel iteration.  40  4.2. An Improved Channel Estimation Algorithm  4.2 An Improved Channel Estimation Algorithm The above discussion is based on the linear approximation (3.20a) near the ML estimation result. Algorithm 1 is, however, a nonlinear estimator. Using the overrelaxation method, we can propose an improved new iterative ML algorithm as shown in Table 4.1. Note that the function L(hl , fl ) in Step 1 is the cost function defined in (3.19) and η = 4/3, which is the asymptotically optimal value given by (4.17). Our numerical results will show that the improved algorithm can achieve good performance, both in terms of MSE and convergence rate. Table 4.1: Algorithm 2 (0)  Initialization: Initialize hl  (0)  and fl , for l = 1, 2, . . . , L. Set the iteration counter i = 0.  Iterations: Step 1: For each l = 1, 2, . . . , L, calculate (i+1) (i) h˜ l = arg min L hl , fl hl  (i+1) (i+1) f˜l = arg min L hl , fl fl  Let  (i+1) hl  (i+1) (i+1) (i) (i+1) (i) + (1 − η )hl and fl + (1 − η ) fl . = η f˜l = η h˜ l (i+1)  Step 2: If max l  hl  (i)  hl  (i+1)  (i)  −hl  fl  × 100% > δh or max  (i)  − fl  × 100% > δ f , let  (i)  fl  l  (i+1)  i = i + 1 and go to Step 1, otherwise output hl  (i+1)  ’s and fl  ’s.  4.3 Numerical Results and Discussions In this section, we present some numerical results by considering the same discrete baseband OFDM uplink transmission discussed in Chapter 3 in our simulation. All the simulation parameters adopted here are the same as those considered in Section 3.5. Here we compare the MSE performance and convergence performance of Algorithm 1 and Algorithm 2 to demonstrate the performance of the improved algorithm. Fig. 4.1 plots the average MSE performance of hˆ l for the two proposed algorithms and Figs. 4.24.4 plot MSE performance of fˆl for the same system with three different fl values. The MSEs are obtained by averaging over 2000 trials, and the theoretical performance predictions of hˆ l and fˆl are 41  4.3. Numerical Results and Discussions obtained using (3.36) and (3.37), respectively. From Figs. 4.1-4.4, we observe that the simulated MSE performance of hˆ l and fˆl for Algorithm 1 and Algorithm 2, along with their theoretical predictions, have excellent agreement over a wide range of SNR values. We can notice a slight difference between the MSE performance of Algorithm 1 and Algorithm 2 when the SNR values are greater than 30 dB. However, the MSEs of both algorithms in this SNR region are less than −50 dB and the differences are therefore negligible. We can conclude that the two algorithms have almost the same performance in terms of MSE. Fig. 4.5 plots the average number of iterations required by the two proposed algorithms to achieve the MSE performance shown in Figs. 4.1-4.4. As expected, Fig. 4.5 shows that the average number of iterations of both algorithms decreases with increasing value of SNR up to 20 dB. The average number of iterations of both algorithms is unchanged when the average SNR value is greater than 20 dB. It is shown in Fig. 4.5 that Algorithm 2 converges faster than Algorithm 1, and it achieves 63% improvement in the number of iterations in the large SNR region. From Fig. 4.5, it is seen that the asymptotically optimal η value obtained under a high SNR assumption also works when the SNR value is low. Figs. 4.6 and 4.7 respectively plot the average MSE performance of hˆ l and fˆl versus the number of iterations. For both figures, the average SNR value is fixed at 25 dB, and we observe that the MSE for Algorithm 2 decreases faster than that of Algorithm 1. Our numerical results clearly demonstrate that Algorithm 2, which uses an SOR technique, can achieve faster convergence than Algorithm 1.  42  4.3. Numerical Results and Discussions  0 Algorithm 1 Algorithm 2 Theoretical value  ˆl (dB) Average MSE of h  −10  −20  −30  −40  −50  −60 −10  0  10 20 Average SNR (dB)  30  40  Figure 4.1: An average MSE performance comparison between Algorithm 1 and Algorithm 2, along with the theoretical values for hˆ l , when δh = 1% and δ f = 1%.  43  4.3. Numerical Results and Discussions  −10 Algorithm 1 Algorithm 2 Theoretical value  −20  MSE of fˆl (dB)  −30  −40  −50  −60  −70 −10  0  10 20 Average SNR (dB)  30  40  Figure 4.2: An MSE performance comparison between Algorithm 1 and Algorithm 2, along with the theoretical values for fˆl , when fl = 0.02, δh = 1% and δ f = 1%.  44  4.3. Numerical Results and Discussions  −10 Algorithm 1 Algorithm 2 Theoretical value  −20  MSE of fˆl (dB)  −30  −40  −50  −60  −70 −10  0  10 20 Average SNR (dB)  30  40  Figure 4.3: An MSE performance comparison between Algorithm 1 and Algorithm 2, along with the theoretical values for fˆl , when fl = 0.04, δh = 1% and δ f = 1%.  45  4.3. Numerical Results and Discussions  −10 Algorithm 1 Algorithm 2 Theoretical value  −20  MSE of fˆl (dB)  −30  −40  −50  −60  −70 −10  0  10 20 Average SNR (dB)  30  40  Figure 4.4: An MSE performance comparison between Algorithm 1 and Algorithm 2, along with the theoretical values for fˆl , when fl = 0.06, δh = 1% and δ f = 1%.  46  4.3. Numerical Results and Discussions  16 Algorithm 1 Algorithm 2  Average Number of Iterations  15 14 13 12 11 10 9 8  0  5  10  15 20 25 Average SNR (dB)  30  35  40  Figure 4.5: An average number of iterations comparison between Algorithm 1 and Algorithm 2 when η = 4/3, δh = 1% and δ f = 1%.  47  4.3. Numerical Results and Discussions  5 Algorithm 1 Algorithm 2  0  ˆl (dB) Average MSE of h  −5 −10 −15 −20 −25 −30 −35 −40  0  5  10 Number of Iterations  15  20  Figure 4.6: The average MSE performance versus number of iterations for hˆ l of Algorithm 1 and Algorithm 2 when η = 4/3 and SNR = 25.  48  4.3. Numerical Results and Discussions  −20 Algorithm 1 Algorithm 2  Average MSE of fˆl (dB)  −25  −30  −35  −40  −45  −50  −55  0  5  10 Number of Iterations  15  20  Figure 4.7: The average MSE performance versus number of iterations for fˆl of Algorithm 1 and Algorithm 2 when η = 4/3 and SNR = 25 dB.  49  Chapter 5 Conclusions In this chapter, we first summarize the contributions of this work, and then suggest some future work.  5.1 Summary of Contributions Through the numerical simulations and mathematical analyses, the performance of our proposed algorithms have been verified. In the simulation, we do not compare the performance of our proposed algorithms with any other solutions because of the uniqueness of our channel model. However, we do compare the simulated estimate’s performance with the CRB, which is the best one can achieve. The contributions of this thesis can be summarized as follows. 1. A discrete CIR model, which allows the CIR to vary within one OFDM symbol, has been developed for OFDM uplink transmission in macrocellular systems. The CIR is characterized by the channel parameters hl and fl , where hl is the complex amplitude for the lth multipath and fl is the corresponding Doppler frequency shift normalized to the duration of an OFDM symbol. In this case, the channel state can be determined by estimating the unknown channel parameters so that the channel estimation problem can be simplified. 2. Taking into account the ICI component, we formulate the channel estimation problem based on our proposed discrete CIR model. Due to the fact that fl is typically less than 0.1 in practical mobile communication scenarios, the ICI factor R(k , k, fl ) is represented by a truncated Taylor series expansion to facilitate the channel estimation problem. We also demonstrate that a second-order Taylor series approximation is sufficient for a typical channel estimation problem. 50  5.2. Future work 3. We propose two iterative ML-based channel estimation algorithms, which do not require channel statistics, for OFDM systems in dispersive time-varying fading channels. These two algorithms can estimate the desired parameters accurately in order to provide channel state information for the receiver. Our proposed channel estimation algorithms are particularly useful for an OFDM system which has already compensated for the frequency offset due to local oscillator mismatch. The MSE performance of the proposed algorithms in large SNR regions has been analyzed using a small perturbation technique. Since our algorithms are iterative, the convergence performance is also analyzed. Based on the analysis, an SOR method is adopted to accelerate the proposed algorithm. The performance analyses are verified by our simulated results.  5.2 Future work In this work, we have proposed two iterative ML-based channel estimation algorithms and demonstrated their performance in terms of MSE. However, channel estimation is only one part of the receiver design. Based on accurate channel estimates, we will design a receiver structure based on the proposed channel model. The bit-error rate (BER) performance of such a system will also be investigated. Single-carrier frequency division multiple access (SC-FDMA) is an OFDM based technique that has been adopted for the Long Term Evolution (LTE) uplink transmission by the third Generation Partnership Project (3GPP). We believe that our proposed channel model and channel estimation algorithms can be applied to such a system. In our future research, we will study the channel estimation problem for SC-FDMA systems using the dispersive, time-varying channel model proposed in this thesis.  51  Bibliography [1] M. K. Ozdemir and H. Arslan, “Channel estimation for wireless OFDM systems,” IEEE Communications Surveys and Tutorials, vol. 9, no. 2, pp. 18–48, July 2007. [2] I. Koffman and V. Roman, “Broadband wireless access solutions based on OFDM access in IEEE 802.16,” IEEE Communications Magazine, vol. 40, no. 4, pp. 96–103, Apr. 2002. [3] I. Barhumi, G. Leus, and M. Moonen, “Optimal training design for MIMO-OFDM systems in mobile wireless channels,” IEEE Transactions on Signal Processing, vol. 51, no. 6, pp. 1615–1624, June 2003. [4] Y. Li and G. L. St¨uber, Orthogonal Frequency Division Multiplexing for Wireless Communications. New York: Springer Publications, 2006. [5] J. G. Proakis, Digital Communications, 4th ed. New York: McGraw-Hill, 2000. [6] J.-J. van de Beek, O. Edfors, M. Sandell, S. Wilson, and P. B¨orjesson, “On channel estimation in OFDM systems,” in Proc. IEEE Vehicular Technology Conference, VTC 1995, Chicago, IL, USA, July 25–28 1995, vol. 2, pp. 815–819. [7] O. Edfors, M. Sandell, J.-J. van de Beek, S. K. Wilson, and P. O. B¨orjesson, “OFDM channel estimation by singular value decomposition,” IEEE Transactions on Communications, vol. 46, no. 7, pp. 931–939, July 1998. [8] P. Chen and H. Kobayashi, “Maximum likelihood channel estimation and signal detection for OFDM systems,” in Proc. IEEE International Conference on Communications, ICC 2002, New York, NY, USA, May 2002, vol. 3, pp. 1640–1645.  52  Chapter 5. Bibliography [9] Y. Li, L. J. Cimini, and N. R. Sollenberger, “Robust channel estimation for OFDM systems with rapid dispersive fading channels,” IEEE Transactions on Communications, vol. 46, no. 7, pp. 902–915, July 1998. [10] J. K. Moon and S. I. Choi, “Performance of channel estimation methods for OFDM systems in a multipath fading channels,” IEEE Transactions on Consumer Electronics, vol. 46, no. 1, pp. 161–170, Feb. 2000. [11] Y. Li, “Pilot-symbol-aided channel estimation for OFDM in wireless systems,” IEEE Transactions on Vehicular Technology, vol. 49, no. 4, pp. 1207–1215, July 2000. [12] Y. Zhao and A. Huang, “A novel channel estimation method for OFDM mobile communications systems based on pilot signals and transform-domain processing,” in Proc. IEEE Vehicular Technology Conference, VTC 1997, Phoenix, AZ, USA, May 4–7 1997, vol. 3, pp. 2089–2093. [13] M. X. Chang and Y. T. Su, “Model-based channel estimation for OFDM signals in Rayleigh fading,” IEEE Transactions on Communications, vol. 50, no. 4, pp. 540–544, Apr. 2002. [14] F. Z. Merli and G. M. Vitetta, “Iterative ML-based estimation of carrier frequency offset, channel impulse response and data in OFDM transmissions,” IEEE Transactions on Communications, vol. 56, no. 3, pp. 497–506, Mar. 2008. [15] G. L. St¨uber, Principles of Mobile Communications. Norwell, MA: Kluwer, 1996. [16] T. S. Rappaport, Wireless Communications: Principles and Practice, 2nd ed. Upper Saddle River: Prentice-Hall, 2001. [17] A. Goldsmith, Wireless Communications. New York: Cambridge University Press, 2005. [18] D. Tse and P. Viswanath, Fundamentals of Wireless Communication. New York: Cambridge University Press, 2005. [19] A. Paulraj, R. Nabar, and D. Gore, Introduction to Space-Time Wireless Communications. Cambridge, UK: Cambridge Univeristy Press, 2003. 53  Chapter 5. Bibliography [20] R. van Nee and R. Prasad, OFDM for Wireless Multimedia Communications. Boston: Artech House, 2000. [21] R. R. Mosier and R. G. Clabaugh, “Kineplex, a bandwidth efficient binary transmission system,” IEEE Transactions on Power Apparatus and Systems, vol. 76, pp. 723–728, Jan. 1958. [22] R. W. Chang, “Synthesis of band limited orthogonal signals for multichannel data transmission,” Bell System Technical Journal, vol. 45, no. 12, pp. 1775–1796, Dec. 1966. [23] B. R. Salzberg, “Performance of an efficient parallel data transmission system,” IEEE Transactions on Communication Technology, vol. COM-15, no. 6, pp. 805–813, Dec. 1967. [24] S. B. Weinstein and P. M. Ebert, “Data transmission by frequency-division multiplexing using the discrete Fourier transform,” IEEE Transactions on Communication Technology, vol. COM-19, no. 5, pp. 628–634, Oct. 1971. [25] A. Gorokhov and J.-P. Linnartz, “Robust OFDM receivers for dispersive time-varying channels: equalization and channel acquisition,” IEEE Transactions on Communications, vol. 52, no. 4, pp. 527–583, Apr. 2004. [26] R. Steele and L. Hanzo, Mobile Radio Communications. New York: Wiley-IEEE Press, 1999. [27] X. Ma, H. Kobayashi, and S. C. Schwartz, “Joint frequency offset and channel estimation for OFDM,” in Proc. IEEE Global Telecommunications Conference, GLOBECOM 2003, San Francisco, CA, USA, Dec. 1–5 2003, vol. 1, pp. 15–19. [28] T. Cui and C. Tellambura, “Joint frequency offset and channel estimation for OFDM systems using pilot symbols and virtual carriers,” IEEE Transactions on Wireless Communications, vol. 6, no. 4, pp. 1193–1202, Apr. 2007. [29] K. B. Petersen and M. S. Petersen, The Matrix Cookbook. http://matrixcookbook.com, 2008.  54  Chapter 5. Bibliography [30] J. A. Fessler and A. O. Hero, “Space-alternating generalized expectation-maximization algorithm,” IEEE Transactions on Signal Processing, vol. 42, no. 10, pp. 2664–2677, Oct. 1994. [31] K. Li and H. Liu, “Joint channel and carrier offset estimation in CDMA communications,” IEEE Transactions on Signal Processing, vol. 47, no. 7, pp. 1811–1822, July 1999. [32] P. Stoica and O. Besson, “Training sequence design for frequency offset and frequencyselective channel estimation,” IEEE Transactions on Communications, vol. 51, no. 11, pp. 1910–1917, Nov. 2003. [33] P. H. Moose, “A technique for orthogonal frequency division multiplexing offset correction,” IEEE Transactions on Communications, vol. 42, no. 10, pp. 2908–2914, Oct. 1994. [34] T. M. Schmidl and D. C. Cox, “Robust frequency and timing synchronization for OFDM,” IEEE Transactions on Communications, vol. 45, no. 12, pp. 1613–1621, Dec. 1997. [35] D. M. Young, Iterative Solution of Large Linear Systems. New York: Academic Press, 1971. [36] G. T. RAN, “Spatial channel model for Multiple Input Multiple Output (MIMO) simulations,” Technical Report, Dec. 2008. [37] G. Cardano and T. R. Witmer, The Great Art, or: The Rules of Algebra, 1st ed. Cambridge: MIT Press, 1968.  55  Appendix A Derivation of the Discrete CIR Model For a multipath CIR described in (3.2), each resolvable path consists of a number of non-resolvable components so that γl (t) can be expressed as  γl (t) = ∑ αl,i e− jφl,i (t)  (A.1)  i  where αl,i and φl,i (t) are respectively the complex amplitude and the phase shift incurred by the ith non-resolvable component of the lth resolvable path. The phase shift φl,i (t) is due to the time delay and Doppler frequency shift and is given by  φl,i (t) = 2π fc τl − 2π fDl,i t  (A.2)  where fc is the carrier frequency and fDl,i is the Doppler frequency shift at the ith non-resolvable component of the lth resolvable path. In this work, we consider an uplink transmission in a macrocellular system where the angular spread is small for those non-resolvable components [36]. Consequently, we can assume fDl,i = fDl and  φl,i (t) = 2π fc τl − 2π fDl t.  (A.3)  Applying (A.3) to (A.1), one can obtain  γl (t) = e j2π fDl t ∑ αl,i e− j2π fc τl  (A.4)  i hl  56  Appendix A. Derivation of the Discrete CIR Model and rewrite (3.2) as L  h(τ ,t) =  ∑ hl e j2π fDl t δ (τ − τl ).  (A.5)  l=1  Therefore, the discrete-time CIR can be expressed as  h[m, n] = h(mTsa , nTsa ) L  =  ∑ hl e l=1 L  =  ∑ hl e l=1  j2π fD Ts nTsa l NTsa  j2π fl n N  δ (mTsa − Tsa τl /Tsa )  δ [m − nl ]  (A.6)  where Tsa is the sampling interval defined as Tsa = Ts /N and Ts is the OFDM symbol duration, fl is the Doppler frequency shift at the lth multipath normalized by 1/Ts , i.e., fl = fDl Ts , and nl is the corresponding delay in samples, i.e., nl = τl /Tsa .  57  Appendix B Derivation of (3.20) Recalling (3.8) and (3.13), the cost function we try to minimize is L  J(hl , fl ) = Y − ∑ hl l=1 L  = Y − ∑ hl l=1  2  Al + Bl fl + Cl fl2 + · · · H  Al + Bl fl + Cl fl2 + · · ·  L  Y − ∑ hl Al + Bl fl + Cl fl2 + · · · l=1  (B.1)  where (·)H denotes the Hermitian transpose of the argument matrix. According to (30), (31), and (32) in [31], we can obtain    ℜ {∆hi }  ∆θ i =   ℑ {∆hi }  ∆ fi   and      Ri j =       ni = −    ∂ 2 J(hl , fl ) ∂ ℜ{hi }∂ ℜ{h j } ∂ 2 J(hl , fl ) ∂ ℑ{hi }∂ ℜ{h j } ∂ 2 J(hl , fl ) ∂ fi ∂ ℜ{h j }  ∂ J(hl , fl ) ∂ ℜ{hi } ∂ J(hl , fl ) ∂ ℑ{hi } ∂ J(hl , fl ) ∂ fi    (B.2)    (B.3)              ∂ 2 J(hl , fl ) ∂ ℜ{hi }∂ ℑ{h j } ∂ 2 J(hl , fl ) ∂ ℑ{hi }∂ ℑ{h j } ∂ 2 J(hl , fl ) ∂ fi ∂ ℑ{h j }  ∂ 2 J(hl , fl ) ∂ ℜ{hi }∂ f j ∂ 2 J(hl , fl ) ∂ ℑ{hi }∂ f j ∂ 2 J(hl , fl ) ∂ fi ∂ f j         (B.4)  where ℜ{·} and ℑ{·} respectively denote the real part and imaginary part of the argument, and ∆hi = hˆ i − hi and ∆ fi = fˆi − fi are the perturbation of tne estimates caused by Gaussian noise ni .  58  Appendix B. Derivation of (3.20) In the following, we will calculate the elements of ni . Using (B.1), one can obtain  ∂ J(hl , fl ) = − Ai + Bi fi + Ci fi2 + · · · ∂ ℜ{hi } L  + Y − ∑ hl l=1  L  H  Y − ∑ hl Al + Bl fl + Cl fl2 + · · · l=1  H  Al + Bl fl + Cl fl2 + · · ·  − Ai + Bi fi + Ci fi2 + · · ·  . (B.5)  Recalling from (3.8), we have L  w = Y − ∑ hl Al + Bl fl + Cl fl2 + · · ·  (B.6)  l=1  and can rewrite (B.5) as  ∂ J(hl , fl ) = − Ai + Bi fi + Ci fi2 + · · · ∂ ℜ{hi }  H  w + wH − Ai + Bi fi + Ci fi2 + · · ·  .  (B.7)  Using the assumption that fl ’s are small, we can neglect the terms containing fi in (B.7) and further simplify (B.7) as  ∂ J(hl , fl ) ≈ −Ai ∂ ℜ{hi }  H  = −wH Ai  w + wH −Ai H  + −wH Ai  = −2ℜ wH Ai .  (B.8)  In a similar way, we can obtain  ∂ J(hl , fl ) = 2ℑ wH Ai ∂ ℑ{hi }  (B.9)  ∂ J(hl , fl ) = −2ℜ hi wH Bi ∂ fi  (B.10)  and  59  Appendix B. Derivation of (3.20) so that we obtain      ni = −     ∂ J(hl , fl ) ∂ ℜ{hi } ∂ J(hl , fl ) ∂ ℑ{hi } ∂ J(hl , fl ) ∂ fi      wH Ai   ℜ     = 2  −ℑ wH A i     ℜ hi wH Bi      .    (B.11)  Based on the derivation above, it is straightforward to calculate the elements of Ri j in a similar way and obtain    AHj Ai   ℜ  H Ri j = 2   −ℑ A j Ai  ℜ hi AHj Bi  ℑ  AHj Ai  ℜ  AHj Ai  ℑ hi AHj Bi  ℜ −ℑ  h∗j BHj Ai h∗j BHj Ai  ℜ{hi h∗j BHj Bi }         (B.12)  where {·}∗ denotes the complex conjugate of the argument. Applying (B.11) and (B.12) to (32) in [31], we can obtain (3.20).  60  Appendix C Derivation of (3.30c) From (3.7) and (3.8), we have BH l Bl =  N−1  ∑ BHl [k]Bl [k].  (C.1)  k=0  Applying (3.9b) to (C.1) and using the uncorrelation property of white training sequence , one obtains BH l Bl  π 2 (N − 1)2 N−1 4π 2 N−1 N−1 1 2 |X[k]| + X[k ] = ∑ ∑ ∑ 2 2 N N k=0 k =k |1 − ω |2 k=0  2  (C.2)  k =0  where ω = exp  j2π (k −k) N  . Assuming |X[k]| =  √1 , N  for k = 0, 1, . . . , N − 1, we can simplify (C.2)  to be BH l Bl  1 π 2 (N − 1)2 4π 2 N−1 N−1 + 3 ∑ ∑ . = 2 N N k=0 k =k |1 − ω |2  (C.3)  k =0  For the sake of simplicity, we denote N−1  s[k] =  =  1  ∑  2 k =k |1 − ω | k =0 N−1  1  ∑ k =k k =0  1−e  j2π (k −k) N  2  ,  k = 0, 1, . . . , N − 1  (C.4)  and rewrite (C.3) as BH l Bl =  π 2 (N − 1)2 4π 2 N−1 + 3 ∑ s[k]. N2 N k=0  (C.5)  When N is fixed, exploiting the fact that s[k] is constant for k = 0, 1, . . . , N − 1, one can obtain N−1  s[k] = s =  ∑ n=1  1 1−e  j2π n N  2  ,  k = 0, 1, . . . , N − 1  (C.6)  61  Appendix C. Derivation of (3.30c) and further simplify (C.5) to be  BH l Bl =  π 2 (N − 1)2 4π 2 + 2 s. N2 N  (C.7)  Using Parseval’s relationship, one can show N−1  s=  ∑ n=1  1 1−e  j2π n N  2  =  1 N−1 1 ∑ sin2 nπ 4 n=1 N  =  (N − 1)(N + 1) . 12  (C.8)  Finally substituting of (C.8) into (C.7), one obtains (3.30c) as desired, i.e.,  BH l Bl =  2π 2 (N − 1)(2N − 1) . 3N 2  (C.9)  62  Appendix D Derivation of the Absolute Value for the Largest Eigenvalue of Mη−1Nη According to algebraic theory, the roots of a cubic equation  x3 + α x2 + β x + γ = 0  (D.1)  are given by [37] x = v−u−  α 3  (D.2)  where u=  3  q2 p3 + 4 27  q ± 2  v=  p 3u  (D.3b)  α2 3  (D.3c)  2α 3 − 9αβ . 27  (D.3d)  p=β−  q=γ+  (D.3a)  After some manipulations, one can calculate that the roots for the characteristic equation in (4.12) are given by p=−  3(1 − η ) + η 2 kl 2 η kl 3  (D.4a)  63  Appendix D. Derivation of the Absolute Value for the Largest Eigenvalue of Mη−1 Nη and q=−  1 9(1 − η ) + 2η 2 kl η 4 kl2 . 27  (D.4b)  To minimize the largest absolute value of the roots of (4.12), we let  δ=  q2 p3 + =0 4 27  (D.5)  and have three equal roots. After applying (D.4) into (D.5), one can show that η is required to satisfy  or  η 2 kl − 4 η + 4 = 0  (D.6)  √ 2(1 − 1 − kl ) . η= kl  (D.7)  After substitution of η into u and v, it can be shown that  u=  3  2 q = − (1 − η ) 2 3  2 v = (1 − η ) 3 x = v−u−  α = α = 1 − η. 3  (D.8a)  (D.8b)  (D.8c)  Thus, the absolute value for the largest eigenvalue of Mη−1 Nη is √ 2(1 − 1 − kl ) . |λ | = |1 − η | = 1 − kl  (D.9)  64  Appendix E Matlab Code %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%  %%  %%  Main F u n c t i o n  %%  %% %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clc ; clear ; %%  Parameter I n i t i a l i z a t i o n  %%  REPEAT = 2000;  % Number o f r e p e a t times , f o r purpose o f average  N f f t = 256;  % Size o f FFT  L = 3;  % L−path channel  % f l = ones ( 1 , L ) ;  % Normalized Doppler f r e q u e n c y ( 1 by L v e c t o r )  % h l = zeros ( 1 , L ) ;  % Amplitude f o r t h e path ( 1 by L v e c t o r )  nl = 0 : (L − 1);  % Delay o f t h e path ( 1 by L v e c t o r )  n modu type = 1 ;  % The t y p e o f m o d u l a t i o n scheme % 1=BPSK, 2=QPSK, 4=16QAM, 6=64QAM  Tx = 1 ;  % Tx=1 => T r a n s m i t o r ; Tx=0 => Receiver  SNR = −10: 5 : 4 0 ; ( − 1.3362 − 0.6918 j )  h l = [(0.2944+1.6236 j )  (0.7143+0.8580 j ) ] ;  f l = [ 0 . 0 2 0.04 0 . 0 6 ] ; delta h 2 = 0.01;  % Stoping threshold  de l t a f 2 = 0.01;  % Stoping threshold  h l i n i t i a l = ones ( 2 , L ) ; f l i n i t i a l = ones ( 2 , L ) ; h l i n i t i a l = (1 + 1 j ) ∗ h l i n i t i a l ;  % I n i t i a l value  f l i n i t i a l = 0.05 ∗ f l i n i t i a l ;  % I n i t i a l value  %%  For A l g o r i t h m 1  %%  h l h a t a v e p 2 A 1 = zeros ( l e n g t h (SNR) , L ) ; f l h a t a v e p 2 A 1 = zeros ( l e n g t h (SNR) , L ) ; hl MSE p2 A1 = zeros ( l e n g t h (SNR) , L ) ; fl MSE p2 A1 = zeros ( l e n g t h (SNR) , L ) ; h l e r r a v e p 2 A 1 = zeros ( l e n g t h (SNR) , L ) ;  65  Appendix E. Matlab Code f l e r r a v e p 2 A 1 = zeros ( l e n g t h (SNR) , L ) ; i t e r a t i o n p 2 A 1 = zeros ( l e n g t h (SNR) , 1 ) ;  %%  For A l g o r i t h m 2  %%  h l h a t a v e p 2 A 2 = zeros ( l e n g t h (SNR) , L ) ; f l h a t a v e p 2 A 2 = zeros ( l e n g t h (SNR) , L ) ; hl MSE p2 A2 = zeros ( l e n g t h (SNR) , L ) ; fl MSE p2 A2 = zeros ( l e n g t h (SNR) , L ) ; h l e r r a v e p 2 A 2 = zeros ( l e n g t h (SNR) , L ) ; f l e r r a v e p 2 A 2 = zeros ( l e n g t h (SNR) , L ) ; i t e r a t i o n p 2 A 2 = zeros ( l e n g t h (SNR) , 1 ) ;  %%  For t h e o r e t i c a l v a l u e  %%  nvar = zeros ( l e n g t h (SNR) , 1 ) ;  % Noise v a r i a n c e  t h e o r e t i c f l = zeros ( l e n g t h (SNR) , ( L + 1 ) ) ;  f o r s n r i n d e x = 1 : l e n g t h (SNR) %%  For A l g o r i t h m 1  %%  h l h a t t e m p p 2 A 1 = zeros ( 1 , L ) ; f l h a t t e m p p 2 A 1 = zeros ( 1 , L ) ; hl MSE temp p2 A1 = zeros ( 1 , L ) ; fl MSE temp p2 A1 = zeros ( 1 , L ) ; h l e r r t e m p p 2 A 1 = zeros ( 1 , L ) ; f l e r r t e m p p 2 A 1 = zeros ( 1 , L ) ; Iteration temp p2 A1 = 0;  %%  For A l g o r i t h m 2  %%  h l h a t t e m p p 2 A 2 = zeros ( 1 , L ) ; f l h a t t e m p p 2 A 2 = zeros ( 1 , L ) ; hl MSE temp p2 A2 = zeros ( 1 , L ) ; fl MSE temp p2 A2 = zeros ( 1 , L ) ; h l e r r t e m p p 2 A 2 = zeros ( 1 , L ) ; f l e r r t e m p p 2 A 2 = zeros ( 1 , L ) ; Iteration temp p2 A2 = 0;  %%  For t h e o r e t i c a l v a l u e  nvar temp = 0 ;  %%  % Noise v a r i a n c e  f o r r e p e a t i n d e x = 1 : REPEAT %%  Generate Data  %%  X b i t = round ( rand ( 1 , ( N f f t ∗ n modu type ) ) ) ; %X b i t = r a n d i n t ( 1 , N f f t ∗ n modu type ) ;  66  Appendix E. Matlab Code %%  Mo d u l a ti o n ( Gray Mapping )  %%  X = mapper ( X b i t , n modu type , Tx ) ; X = X / sqrt ((X ∗ X ’ ) ) ;  %%  %%% Normalize X  Calculate Coefficiences  %%  [ Al , Bl , Cl , Dl ] = c a l c u c o e f f ( N f f t , L , n l , X ) ;  %%  C a l c u l a t e t h e o u t p u t o f FFT (Y [ k ] ) u s i n g Eq . 4  %%  [ Y , n v a r i ] = c a l c u c o r r u p t e d ( N f f t , L , f l , h l , n l , X ,SNR( s n r i n d e x ) ) ; nvar temp = nvar temp + n v a r i ;  %%  Channel E s t i m a t i o n u s i n g A l g o r i t h m 1  %%  temp iter p2 = 0; p = 2;  % The o r d e r o f t h e T a y l o r expansion . % I t can n o t work p r o p e r l y when p =1.  [ hl hat ave p2 A1 ( snr index , : ) , fl hat ave p2 A1 ( snr index , : ) , i t e r a t i o n p 2 A 1 ( snr index ) ] . . . = c h a n n e l e s t A 1 ( N f f t , L , p , Al , Bl , Cl , Dl , Y , d e l t a h 2 , d e l t a f 2 , h l i n i t i a l ,  f l i n i t i a l );  i t e r a t i o n p 2 A 1 ( s n r i n d e x ) = 1000 − i t e r a t i o n p 2 A 1 ( s n r i n d e x ) ; %%  Channel E s t i m a t i o n u s i n g A l g o r i t h m 2  %%  [ hl hat ave p2 A2 ( snr index , : ) , fl hat ave p2 A2 ( snr index , : ) , i t e r a t i o n p 2 A 2 ( snr index ) ] . . . = c h a n n e l e s t A 2 ( N f f t , L , p , Al , Bl , Cl , Dl , Y , d e l t a h 2 , d e l t a f 2 , h l i n i t i a l ,  f l i n i t i a l );  i t e r a t i o n p 2 A 2 ( s n r i n d e x ) = 1000 − i t e r a t i o n p 2 A 2 ( s n r i n d e x ) ; %%  C a l c u l a t e t h e MSE  %%  For A l g o r i t h m 1  %% %%  hl hat temp p2 A1=hl hat temp p2 A1+hl hat ave p2 A1 ( snr index , : ) ; fl hat temp p2 A1=fl hat temp p2 A1+fl hat ave p2 A1 ( snr index , : ) ; hl temp p2 A1 = abs ( ( h l h a t a v e p 2 A 1 ( s n r i n d e x , : ) − h l ) ) . ˆ 2 ; f l t e m p p 2 A 1 = abs ( ( f l h a t a v e p 2 A 1 ( s n r i n d e x , : ) − f l ) ) . ˆ 2 ; hl MSE temp p2 A1 = hl MSE temp p2 A1 + hl temp p2 A1 ; fl MSE temp p2 A1 = fl MSE temp p2 A1 + f l t e m p p 2 A 1 ; Iteration temp p2 A1 = Iteration temp p2 A1 + iteration p2 A1 ( snr index ) ;  %%  For A l g o r i t h m 2  %%  hl hat temp p2 A2=hl hat temp p2 A2+hl hat ave p2 A2 ( snr index , : ) ; fl hat temp p2 A2=fl hat temp p2 A2+fl hat ave p2 A2 ( snr index , : ) ; hl temp p2 A2 = abs ( ( h l h a t a v e p 2 A 2 ( s n r i n d e x , : ) − h l ) ) . ˆ 2 ; f l t e m p p 2 A 2 = abs ( ( f l h a t a v e p 2 A 2 ( s n r i n d e x , : ) − f l ) ) . ˆ 2 ; hl MSE temp p2 A2 = hl MSE temp p2 A2 + hl temp p2 A2 ; fl MSE temp p2 A2 = fl MSE temp p2 A2 + f l t e m p p 2 A 2 ; Iteration temp p2 A2 = Iteration temp p2 A2 + iteration p2 A2 ( snr index ) ; end  67  Appendix E. Matlab Code  %%  For A l g o r i t h m 1  %%  h l h a t a v e p 2 A 1 ( s n r i n d e x , : ) = h l h a t t e m p p 2 A 1 / REPEAT ; f l h a t a v e p 2 A 1 ( s n r i n d e x , : ) = f l h a t t e m p p 2 A 1 / REPEAT ; hl MSE p2 A1 ( s n r i n d e x , : ) = hl MSE temp p2 A1 / REPEAT ; fl MSE p2 A1 ( s n r i n d e x , : ) = fl MSE temp p2 A1 / REPEAT ; hl err ave p2 A1 ( snr index , : ) = hl hat ave p2 A1 ( snr index , : ) − hl ; fl err ave p2 A1 ( snr index , : ) = fl hat ave p2 A1 ( snr index , : ) − f l ; i t e r a t i o n p 2 A 1 ( s n r i n d e x ) = I t e r a t i o n t e m p p 2 A 1 / REPEAT ;  %%  For A l g o r i t h m 2  %%  h l h a t a v e p 2 A 2 ( s n r i n d e x , : ) = h l h a t t e m p p 2 A 2 / REPEAT ; f l h a t a v e p 2 A 2 ( s n r i n d e x , : ) = f l h a t t e m p p 2 A 2 / REPEAT ; hl MSE p2 A2 ( s n r i n d e x , : ) = hl MSE temp p2 A2 / REPEAT ; fl MSE p2 A2 ( s n r i n d e x , : ) = fl MSE temp p2 A2 / REPEAT ; hl err ave p2 A2 ( snr index , : ) = hl hat ave p2 A2 ( snr index , : ) − hl ; fl err ave p2 A2 ( snr index , : ) = fl hat ave p2 A2 ( snr index , : ) − f l ; i t e r a t i o n p 2 A 2 ( s n r i n d e x ) = I t e r a t i o n t e m p p 2 A 2 / REPEAT ;  %%  For t h e o r e t i c a l v a l u e  %%  nvar ( s n r i n d e x ) = nvar temp / REPEAT ; end  %%  Plot Figure  %%  hl temp p2 A1 = zeros ( l e n g t h (SNR) , 1 ) ; f l t e m p p 2 A 1 = zeros ( l e n g t h (SNR) , 1 ) ; t h e o r e t i c h l = nvar ∗ 5 / 2 ; h l n o r m a l = abs ( h l ) . ˆ 2 ; for i = 1 : L t h e o r e t i c f l ( : , i ) = nvar / h l n o r m a l ( i ) ; t h e o r e t i c f l ( : , ( L + 1)) = t h e o r e t i c f l ( : , ( L + 1))+ t h e o r e t i c f l ( : , i ) ; end t h e o r e t i c f l = t h e o r e t i c f l / 2 / p i ˆ2 ∗ 3; t h e o r e t i c f l (: , (L + 1)) = t h e o r e t i c f l (: , (L + 1)) / L; hl temp p2 A2 = zeros ( l e n g t h (SNR) , 1 ) ; f l t e m p p 2 A 2 = zeros ( l e n g t h (SNR) , 1 ) ; for i = 1 : L hl temp p2 A1 = hl MSE p2 A1 ( : , i ) + hl temp p2 A1 ; hl temp p2 A2 = hl MSE p2 A2 ( : , i ) + hl temp p2 A2 ; end hl ave MSE p2 A1 = hl temp p2 A1 ’ / L ; hl ave MSE p2 A2 = hl temp p2 A2 ’ / L ; hl ave MSE p2 A1 = 10 ∗ log10 ( hl ave MSE p2 A1 ) ;  68  Appendix E. Matlab Code fl ave MSE p2 A1 = 10 ∗ log10 ( fl MSE p2 A1 ) ; hl ave MSE p2 A2 = 10 ∗ log10 ( hl ave MSE p2 A2 ) ; fl ave MSE p2 A2 = 10 ∗ log10 ( fl MSE p2 A2 ) ; t h e o r e t i c h l f i n a l = 10 ∗ log10 ( t h e o r e t i c h l ) ; t h e o r e t i c f l f i n a l = 10 ∗ log10 ( t h e o r e t i c f l ) ; for p = 2 : 3 switch p case 2 p l o t (SNR, hl ave MSE p2 A1 , ’ − ko ’ ) ; h o l d on p l o t (SNR, hl ave MSE p2 A2 , ’ −mp ’ ) ; h o l d on p l o t (SNR, t h e o r e t i c h l f i n a l , ’ − bd ’ ) ; g r i d on ; legend ( ’ A l g o r i t h m 1 ’ , ’ A l g o r i t h m 2 ’ , ’ T h e o r e t i c a l value ’ ) ; x l a b e l ( ’ Average SNR ( dB ) ’ ,  ’ interpreter ’ ,  y l a b e l ( ’ Average MSE o f $ \ h a t h l $ ( dB ) ’ ,  ’ latex ’ ) ; ’ interpreter ’ ,  ’ latex ’ ) ;  end end  for i = 1 : L figure p l o t (SNR, fl ave MSE p2 A1 ( : , i ) , ’ − ko ’ ) ; h o l d on p l o t (SNR, fl ave MSE p2 A2 ( : , i ) , ’ −mp ’ ) ; h o l d on p l o t (SNR, t h e o r e t i c f l f i n a l ( : , i ) , ’ − bd ’ ) ; g r i d on ; legend ( ’ A l g o r i t h m 1 ’ , ’ A l g o r i t h m 2 ’ , ’ T h e o r e t i c a l value ’ ) ; x l a b e l ( ’ Average SNR ( dB ) ’ ,  ’ interpreter ’ ,  y l a b e l ( ’MSE o f $ \ h a t f l $ ( dB ) ’ ,  ’ latex ’ ) ;  ’ interpreter ’ ,  ’ latex ’ ) ;  end x d = SNR( 3 : end ) ; y d A1 = i t e r a t i o n p 2 A 1 ( 3 : end ) ; y d A2 = i t e r a t i o n p 2 A 2 ( 3 : end ) ; figure p l o t ( x d , y d A1 , ’ − ko ’ ) ; h o l d on p l o t ( x d , y d A2 , ’ −mp ’ ) ; h o l d on g r i d on legend ( ’ A l g o r i t h m 1 ’ , ’ A l g o r i t h m 2 ’ ) ; x l a b e l ( ’ Average SNR ( dB ) ’ ,  ’ interpreter ’ ,  ’ latex ’ ) ;  69  Appendix E. Matlab Code y l a b e l ( ’ Average Number o f I t e r a t i o n s ’ ,  %  The End  ’ interpreter ’ ,  ’ latex ’ ) ;  %  70  Appendix E. Matlab Code f u n c t i o n data mapped = mapper ( d a t a i n p u t , n modu type , Tx )  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%  %%  %%  F i l e : mapear .m  %%  %%  %%  %%  F u n c t i o n : Map t h e i n p u t b i t s i n t o t h e modulated symbols  %%  %%  a c c o r d i n g t o t h e m o du l a t i o n t y p e decided by  %%  %%  t h e parameter ’ n modu type ’ a t t h e t r a n s m i t o r .  %%  %%  Or demap t h e r e c e i v e d symbols i n t o b i t s a t t h e  %%  %%  receiver .  %%  %%  %% Parameters : d a t a i n p u t => i n p u t b i t s  %%  %%  %%  n modu type => mo d u la t i o n t y p e  %%  %%  (1=BPSK, 2=QPSK, 4=16QAM, 6=64QAM)  %%  %%  Tx => f l a g i n d i c a t e s t h e t r a n s m i t o r o r r e c e i v e r  %%  %%  (1= t r a n s m i t o r , 0= r e c e i v e r )  %%  %%  %%  %%  Outputs : A m a t r i x o f t h e mapped symbols a c c o r d i n g t o t h e  %%  %%  m o d u l a t i o n scheme used a t t h e t r a n s m i t o r s i d e .  %%  %%  Or a v e c t o r o f t h e demapped b i t s a c c o r d i n g t o t h e  %%  %%  m o d u l a t i o n scheme .  %%  %%  %%  %%  N o t i c e : The l e n g t h o f ’ d a t a i n p u t ’ must be a i n t e g e r m u l t i p l e o f ’ n modu type ’ a t t h e t r a n s m i t o r s i d e .  %% %%  %% %% %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  %%  Decide t h e c o n s t e l l a t i o n f i r s t  %%  [M, M1, M2, type mapping , c ] = p a r a m e t e r s c o n s t e l l a t i o n ( n modu type ) ; a l p h a b e t = b i t s y m b o l (M, type mapping , M1, M2 ) ; i f n modu type ˜= 1 c o n s t e l l a t i o n g r a y = alphabet ( : , 3 ) + 1 j ∗ alphabet ( : , 2 ) ; else c o n s t e l l a t i o n g r a y = [ −1 1 ] ’ ; end l = length ( data input ) ;  %%  R e a l i z e t h e mapping o r demapping  %%  i f Tx==1 m a t r i x d a t a = reshape ( d a t a i n p u t , n modu type , l / n modu type ) ; m data decimal = bi2de ( m a t r i x d a t a ’ ,  ’ l e f t −msb ’ ) ;  f o r i = 1 : ( l / n modu type )  71  Appendix E. Matlab Code v d a t a d e c i m a l = m data decimal ( i ) ; v encoded = genqammod ( v d a t a d e c i m a l , c o n s t e l l a t i o n g r a y ) ; o u t p u t ( : , i ) = v encoded ; end data mapped = c . ∗ o u t p u t ;  e l s e i f Tx==0 data normalized = data input . / c ; for i = 1 : l v data mapped = d a t a n o r m a l i z e d ( i ) ; v decoded = genqamdemod ( v data mapped , c o n s t e l l a t i o n g r a y ) ; d a t a d e c i m a l ( : , i ) = v decoded ; data mapped = de2bi ( d a t a d e c i m a l , n modu type ,  ’ l e f t −msb ’ ) ’ ;  data mapped = data mapped ( : ) ’ ; end end  72  Appendix E. Matlab Code f u n c t i o n [M, M1, M2, type mapping , c ] = p a r a m e t e r s c o n s t e l l a t i o n ( n modu type )  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% %%  %% F i l e : p a r a m e t e r s c o n s t e l l a t i o n .m  %% %%  %% F u n c t i o n : Decide t h e c o n s t e l l a t i o n parameters a c c o r d i n g  %%  t o m o d u l a t i o n scheme t o be used .  %% %%  %% %% %%  Parameters : n modu type => t y p e o f m o du l a t i o n scheme  %%  (1=BPSK, 2=QPSK, 4=16QAM, 6=64QAM)  %% %%  %%  %% %% %%  Outputs : The parameters which are used t o c a l l a n o t h e r  %%  %%  f u n c t i o n to c a l c u l a t e the c o n s t e l l a t i o n using  %%  %%  Gray l a b e l i n g .  %%  %%  %%  %%  N o t i c e : Here we should pay a t t e n t i o n t o t h e parameter ’ c ’  %%  %%  We use i t t o n o r m a l i z e t h e c o n s t e l l a t i o n .  %%  %%  %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  s w i t c h n modu type case 1 type mapping = ’MPSK’ ;  % For BPSK ( As a PAM) % Only t a k e i n t o account A i c  M = 2; M1 = 0 ; M2 = 0 ; c = 1; case 2  % For QPSK ( L i k e 4−QAM)  type mapping = ’QAM’ ; c = 1/ sqrt ( 2 ) ; case 4  % For 16−QAM  type mapping = ’QAM’ ; c = 1/ sqrt (10); case 6  % For 64−QAM  type mapping = ’QAM’ ; c = 1/ sqrt (42); end i f n modu type ˜= 1 M = 0; M1 = s q r t ( 2 ˆ n modu type ) ; M2 = s q r t ( 2 ˆ n modu type ) ; end  73  Appendix E. Matlab Code f u n c t i o n [ a l p h a b e t ] = b i t s y m b o l ( M, type , M1, M2)  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% %%  %% F i l e : b i t s y m b o l .m  %%  %% %%  %% F u n c t i o n : Return t h e a l p h a b e t o f t h e c o n s t e l l a t i o n a c c o r d i n g  %%  %%  t o t h e m o du l a t i on scheme used . Here t h e encoding  %%  %%  i s done by Gray l a b e l i n g  %%  %%  %%  %%  Parameters : M => t o t a l number o f p o i n t s i n t h e c o n s t e l l a t i o n  %%  %%  I t should equal t o 2 ˆ k , where k i s t h e number o f  %%  %%  b i t s t h a t are grouped t o g e t h e r t o form a symbol ,  %%  %%  except QAM, which can come as t h e p r o d u c t o f t h e  %%  %%  arguments M1 and M2.  %%  %%  t y p e => t y p e o f m od u l a t io n scheme  %%  %%  PAM, MPSK, QAM( r e c t a n g u l a r ) .  %%  %%  M1 & M2 => number o f p o i n t s on each a x l e  %%  %% %%  %% Outputs : t h e a l p h a b e t o f t h e c o n s t e l l a t i o n a c c o r d i n g  %%  t o t h e m o d u l a t i o n scheme used .  %% %%  %% %% %%  N o t i c e : When u s i n g QAM, M i s t h e p r o d u c t o f M1 by M2.  %%  Here t h e encoding i s done by Gray l a b e l i n  %%  %% %% %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  %%  Block i n i t i a l i z a t i o n  i f strcmp ( type ,  %%  ’QAM’ )  k1 = c e i l ( l o g 2 (M1) ) ; k2 = c e i l ( l o g 2 (M2) ) ; M1 = 2 ˆ k1 ; M2 = 2 ˆ k2 ; M = M1 ∗ M2; Aicd = zeros ( 1 , k1 ) ; Aisd = zeros ( 1 , k2 ) ; t a b l e 1 = zeros (M1, 2 ) ; t a b l e 2 = zeros (M2, 2 ) ; a l p h a b e t = zeros (M, 3 ) ; d1 = 0 : 1 : M1−1; d1 = d1 ’ ; d2 = 0 : 1 : M2−1; d2 = d2 ’ ;  74  Appendix E. Matlab Code i n d 1 = bi2de ( f l i p l r ( g r a y 2 b i ( f l i p l r ( de2bi ( d1 ) ) ) ) ) ; t a b l e 1 = [ d1 , i n d 1 + 1 ] ; i n d 2 = bi2de ( f l i p l r ( g r a y 2 b i ( f l i p l r ( de2bi ( d2 ) ) ) ) ) ; t a b l e 2 = [ d2 , i n d 2 + 1 ] ; else k = c e i l ( l o g 2 (M) ) ; M = 2 ˆ k; Aicd = zeros ( 1 , k ) ; Aisd = zeros ( 1 , k ) ; t a b l e = zeros (M, 2 ) ; a l p h a b e t = zeros (M, 3 ) ; d = 0 : 1 : M−1; d = d’; i n d = bi2de ( f l i p l r ( g r a y 2 b i ( f l i p l r ( de2bi ( d ) ) ) ) ) ; table = [ d , ind +1]; end  %%  Block computation  i f strcmp ( type ,  %%  ’PAM’ )  Aicd = −(M−1) : 2 : M−1; Aisd = [ ] ; % Use t h e a l p h a b e t f o r i =1 : M i n d e x = f i n d i n d e x ( i −1, t a b l e ) ; a l p h a b e t ( i , : ) = [ i −1, Aicd ( i n d e x ) , 0 ] ; end e l s e i f strcmp ( type ,  ’MPSK’ )  angle = 0:2 ∗ p i /M:2 ∗ p i ∗ (M− 1)/M; Aicd = cos ( angle ) ; Aisd = s i n ( angle ) ; % Use t h e a l p h a b e t f o r i =1 : M i n d e x = f i n d i n d e x ( i −1, t a b l e ) ; a l p h a b e t ( i , : ) = [ i −1, Aicd ( i n d e x ) , Aisd ( i n d e x ) ] ; end e l s e i f strcmp ( type ,  ’QAM’ )  Aicd = −(M1−1) : 2 : M1−1; Aisd = (M2−1) : −2 : −(M2− 1); % Use t h e a l p h a b e t f o r i =1 : M1 f o r j =1 : M2 index1 = f i n d i n d e x ( i −1, t a b l e 1 ) ; index2 = f i n d i n d e x ( j −1, t a b l e 2 ) ;  75  Appendix E. Matlab Code l = i + M1 ∗ ( j − 1); a l p h a b e t ( l , : ) = [ l −1, Aicd ( index1 ) , Aisd ( index2 ) ] ; end end end  76  Appendix E. Matlab Code function b = gray2bi ( g )  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%  %%  %%  GRAY2BI c o n v e r t s Gray encoded sequence i n t o t h e b i n a r y  %%  %%  sequence . I t i s assumed t h a t t h e most s i g n i f i c a n t b i t i s  %%  %%  t h e l e f t hand s i d e b i t .  %%  %%  %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  % copy t h e msb : b(: ,1) = g(: ,1);  for i = 2: size (g , 2 ) , b ( : , i ) = x o r ( b ( : , i − 1) , g ( : , i ) ) ; end  return ;  77  Appendix E. Matlab Code f u n c t i o n i n d i c e = f i n d i n d e x ( num , t a b l e )  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%  %%  %%  F i l e : f i n d i n d e x .m  %%  %% %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  dim = s i z e ( t a b l e ) ; num col = dim ( 1 ) ;  f o r i =1: num col i f num== t a b l e ( i , 1 ) ; indice = table ( i , 2 ) ; return end end  78  Appendix E. Matlab Code f u n c t i o n [ Al , Bl , Cl , Dl ] = c a l c u c o e f f ( N f f t , L , n l , X )  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%  %%  %%  F i l e : c a l c u c o e f f .m  %%  %% %%  %%  F u n c t i o n : C a l c u l a t e t h e c o e f f i c i e n t s o f Equation 8 .  %%  Then o u t p u t A l [ k ] , B l [ k ] , Cl [ k ] and Dl [ k ] .  %%  %% %% %%  %%  Parameters : N f f t => Size o f FFT  %%  %%  L => L−path channel  %%  %%  n l => Delay o f t h e path  %%  %%  X => T r a i n i n g sequence  %%  %%  %%  %%  Outputs : Al , Bl , Cl and Dl  %%  %% %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  A l = zeros ( L , N f f t ) ; for l = 1 : L for k = 1 : Nfft A l ( l , k ) = X ( k ) ∗ exp(−1 j ∗ 2 ∗ p i ∗ ( k − 1 ) ∗ n l ( l ) / N f f t ) ; w = exp ( 1 j ∗ 2 ∗ p i ∗ ( ( 1 : N f f t ) − k + 1e −13) / N f f t ) ; temp = X . ∗ exp(−1 j ∗ 2 ∗ p i ∗ ( 0 : ( N f f t − 1 ) ) ∗ n l ( l ) / N f f t ) . . . . ∗ ( −2 ∗ 1 j ∗ p i / N f f t . / ( 1 − w ) ) ; temp ( k ) = A l ( l , k ) ∗ 1 j ∗ p i ∗ ( N f f t − 1 ) / N f f t ; B l ( l , k ) = sum ( temp ) ;  temp = X . ∗ exp(−1 j ∗ 2 ∗ p i ∗ ( 0 : ( N f f t − 1 ) ) ∗ n l ( l ) / N f f t ) . . . . ∗ ( 2 ∗ p i ˆ 2 . ∗ ( 1 − 2 ∗ w . / (w−1) / N f f t ) . / (1−w) / N f f t ) ; temp ( k ) = −A l ( l , k ) ∗ p i ˆ 2 ∗ ( 2 / 3 − 2 / 3 / N f f t / N f f t − ( N f f t −1) / N f f t / N f f t ) ; Cl ( l , k ) = sum ( temp ) ;  temp = X . ∗ exp(−1 j ∗ 2 ∗ p i ∗ ( 0 : ( N f f t − 1 ) ) ∗ n l ( l ) / N f f t ) . . . . ∗ ( 1 j ∗ 4 ∗ p i ˆ 3 ∗ ( N f f t ˆ 2 ∗ (1−w ) . ˆ 2 + 3 ∗ w . ∗ ( N f f t + 1 . . . + w + N f f t ∗ w ) ) / 3 / N f f t ˆ 3 / (1−w ) . ˆ 3 ) ; temp ( k ) = −A l ( l , k ) ∗ 1 j ∗ p i ˆ 3 ∗ ( N f f t − 1)ˆ2 / 3 / N f f t ˆ 2 ; Dl ( l , k ) = sum ( temp ) ; end end  79  Appendix E. Matlab Code f u n c t i o n [ Y , nvar ] = c a l c u c o r r u p t e d ( N f f t , L , f l , h l , n l , X , SNR)  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%  %%  %%  F i l e : c a l c u c o r r u p t e d .m  %%  %%  %%  %%  F u n c t i o n : C a l c u l a t e t h e c o r r u p t e d s i g n a l , undergoing f a d i n g ,  %%  adding AWGN a c c o r d i n g t o Equation ( 4 ) .  %% %%  %%  %%  %%  Parameters : N f f t => Size o f FFT  %%  %%  L => L−path channel  %%  %%  f l => Normalized Doppler f r e q u e n c y  %%  %%  h l => Complex valued channel g a i n  %%  %%  n l => Delay o f t h e path  %%  %%  X => T r a i n i n g sequence  %%  %%  SNR => SNR v a l u e  %%  %%  %%  %%  Outputs : The c o r r u p t e d s i g n a l a t t h e r e c e i v e r s i d e .  %%  %%  %%  %%  N o t i c e : We do n o t add AWGN i f we j u s t want t h e c o r r u p t e d  %%  %%  s i g n a l f o r c a l c u l a t i n g t h e e q u i v a l e n t SNR t o  %%  %%  q u a n t i f y t h e a p p r o x i m a t i o n . Then we o n l y have  %%  %%  6 i n p u t arguments w i t h o u t SNR.  %%  %%  %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  switch nargin case 6 Y = zeros ( 1 , N f f t ) ; R = zeros ( L , 2∗ N f f t − 1); kk = ( 1 − N f f t ) : ( N f f t − 1 ) ; % C a l c u l a t e t h e r e c e i v e d s i g n a l w i t h o u t AWGN for l = 1 : L R( l , : ) = ( 1 − exp ( j ∗ 2 ∗ p i ∗ f l ( l ) ) )  . / ( 1 − exp ( j ∗ 2 . . .  ∗ p i . ∗ ( kk + f l ( l ) + 1e −13) / N f f t ) ) / N f f t ; end for k = 1 : Nfft for l = 1 : L for k  = 1 : Nfft  Y ( k ) = Y ( k ) + exp(− j ∗ 2 ∗ p i ∗ ( k −1) ∗ n l ( l ) /  ...  N f f t ) ∗ X ( k ) ∗ h l ( l ) ∗ R( l , ( N f f t −k+ k ) ) ; end  80  Appendix E. Matlab Code end end nvar = 0 ; case 7 SNR l = 10 ˆ (SNR / 1 0 ) ; Y = zeros ( 1 , N f f t ) ; R = zeros ( L , 2∗ N f f t − 1); kk = ( 1 − N f f t ) : ( N f f t − 1 ) ; % C a l c u l a t e t h e r e c e i v e d s i g n a l w i t h o u t AWGN for l = 1 : L R( l , : ) = ( 1 − exp ( j ∗ 2 ∗ p i ∗ f l ( l ) ) )  . / ( 1 − exp ( j ∗ 2 . . .  ∗ p i . ∗ ( kk + f l ( l ) + 1e −13) / N f f t ) ) / N f f t ; end for k = 1 : Nfft for l = 1 : L for k  = 1 : Nfft  Y ( k ) = Y ( k ) + exp(− j ∗ 2 ∗ p i ∗ ( k −1) ∗ n l ( l ) /  ...  N f f t ) ∗ X ( k ) ∗ h l ( l ) ∗ R( l , ( N f f t −k+ k ) ) ; end end end  % Add AWGN t o t h e s i g n a l Es = mean ( abs (Y ) . ˆ 2 ) ; nvar = Es / SNR l ; AWGNoise = s q r t ( nvar / 2 ) . ∗ ( randn ( 1 , N f f t ) + 1 j . ∗ randn ( 1 , N f f t ) ) ; Y = Y + AWGNoise ; end  81  Appendix E. Matlab Code f u n c t i o n [ h l h a t , f l h a t , count ] = c h a n n e l e s t A 1 ( N f f t , L , p , Al , Bl , Cl , Dl , Y , d e l t a h , d e l t a f , h l i n i t i a l ,  ...  fl initial )  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%  %%  %% F i l e : c h a n n e l e s t A 1 .m  %%  %%  %%  %% F u n c t i o n : E s t i m a t e t h e channel parameters ( h l & f l ) u s i n g A l g o r i t h m 1.%% %%  %%  %% Parameters :  %%  %%  N f f t => Size o f FFT  %%  %%  L => L−path channel  %%  %%  p => t h e o r d e r o f T a y l o r S e r i e s Expansion  %%  %%  Al , Bl , Cl , Dl => t h e c o e f f i c i e n t s o f t h e a p p r o x i m a t i o n  %%  %%  Y => t h e r e c e i v e d s i g n a l  %%  %%  d e l t a h => t h e t h r e s h o l d o f h l ( s t o p s t h e i t e r a t i o n )  %%  %%  d e l t a f => t h e t h r e s h o l d o f f l  %%  %%  h l i n i t i a l => t h e i n i t i a l v a l u e o f e s t i m a t e ’ h l h a t ’  %%  %%  f l i n i t i a l => t h e i n i t i a l v a l u e o f e s t i m a t e ’ f l h a t ’  %%  %%  %%  %% Outputs : The esimates ’ h l h a t ’ ,  ’ f l h a t ’ and t h e i t e r a t i o n  ’ count ’  %%  %% %%  %% N o t i c e : The o r d e r o f accuracy depends on N f f t and p . We have  %%  %%  d i f f e r e n t accuracy f o r d i f f e r e n t path . The i t e r a t i o n t i m e  %%  %%  depends on t h e i n i t i a l v a l u e and t h e t h r e s h o l d .  %%  %%  %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  %%% I n i t i a l i z e h l and f l %%% i f n a r g i n == 10 YYY = Y . ’ ; AAA = A l . ’ ; BBB = B l . ’ ; temp H = i n v ( A l ∗ AAA) ∗ A l ∗ YYY ; temp F = i n v ( B l ∗ BBB) ∗ B l ∗ (YYY − AAA ∗ temp H ) . / temp H ; h l i n i t i a l = [ temp H temp H ] ; f l i n i t i a l = [ temp F temp F ] ; hl initial = hl initial . ’; fl initial = fl initial . ’; end count = 0 ; h l h a t = zeros ( 1 , L ) ; f l h a t = zeros ( 1 , L ) ;  82  Appendix E. Matlab Code switch p case 2  % second o r d e r a p p r o x i m a t i o n  i t e r a t i o n = 1000; while ( i t e r a t i o n ) for l = 1 : L Alpha h = zeros ( 1 , N f f t ) ; A l p h a f = zeros ( 1 , N f f t ) ; if  l == 1 for k = 2 : L Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k ,  :)...  + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) ) ; end e l s e i f l == L for k = 1 : L − 1 Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k ,  :)...  + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) ) ; end else for k = 1 : l − 1 Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k ,  :)...  + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) ) ; end for k = l + 1 : L Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗  ( Al ( k ,  :)...  + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) ) ; end end % update h l Alpha h = Y − Alpha h ; Beta h = A l ( l , : ) + f l i n i t i a l ( 2 , l ) . ∗ B l ( l , : ) + f l i n i t i a l ( 2 , l ) ˆ 2 . ∗ Cl ( l , : ) ; temp = h l i n i t i a l ( 2 , l ) ; h l i n i t i a l ( 2 , l ) = Beta h ∗ Alpha h ’ / ( Beta h ∗ Beta h ’ ) ; h l i n i t i a l (2 , l ) = conj ( h l i n i t i a l (2 , l ) ) ; h l i n i t i a l ( 1 , l ) = temp ; % update f l A l p h a f = Alpha h − h l i n i t i a l ( 2 , l ) . ∗ A l ( l , : ) ; Beta f = h l i n i t i a l (2 , l ) .∗ Bl ( l , : ) ; E t a f = h l i n i t i a l ( 2 , l ) . ∗ Cl ( l , : ) ; a = Eta f ∗ Eta f ’ ∗ 4; b = ( Beta f ∗ Eta f ’ + Eta f ∗ Beta f ’ ) ∗ 3; c = ( Beta f ∗ Beta f ’ − Alpha f ∗ Eta f ’ − Eta f ∗ Alpha f ’ ) ∗ 2; d = −( A l p h a f ∗ B e t a f ’ + B e t a f ∗ A l p h a f ’ ) ; ppp = [ a b c d ] ;  83  Appendix E. Matlab Code r r r = r o o t s ( ppp ) ; temp = f l i n i t i a l ( 2 , l ) ; f l i n i t i a l (2 , l ) = r r r ( length ( r r r ) ) ; f l i n i t i a l ( 1 , l ) = temp ; end i f ( max ( abs ( h l i n i t i a l ( 2 , : ) − h l i n i t i a l ( 1 , : ) )  . / abs ( h l i n i t i a l ( 1 , : ) ) ) < d e l t a h ) . . .  && ( max ( abs ( f l i n i t i a l ( 2 , : ) − f l i n i t i a l ( 1 , : ) )  . / abs ( f l i n i t i a l ( 1 , : ) ) ) < d e l t a f )  count = i t e r a t i o n ; i t e r a t i o n = 0; else i t e r a t i o n = i t e r a t i o n − 1; end end hl hat = h l i n i t i a l (2 , : ) ; f l h a t = f l i n i t i a l (2 , : ) ;  case 3  % t h i r d order approximation  i t e r a t i o n = 1000; while ( i t e r a t i o n ) for l = 1 : L Alpha h = zeros ( 1 , N f f t ) ; A l p h a f = zeros ( 1 , N f f t ) ; if  l == 1 for k = 2 : L Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k , : ) + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) . . . + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 3 . ∗ Dl ( k , : ) ) ; end  e l s e i f l == L for k = 1 : L − 1 Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k , : ) + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) . . . + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 3 . ∗ Dl ( k , : ) ) ; end else for k = 1 : l − 1 Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k , : ) + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) . . . + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 3 . ∗ Dl ( k , : ) ) ; end for k = l + 1 : L Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k , : ) + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) . . . + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 3 . ∗ Dl ( k , : ) ) ; end end % update h l  84  Appendix E. Matlab Code Alpha h = Y − Alpha h ; Beta h = A l ( l , : ) + f l i n i t i a l ( 2 , l ) . ∗ B l ( l , : ) + f l i n i t i a l ( 2 , l ) ˆ 2 . ∗ Cl ( l ,  :)...  + f l i n i t i a l ( 2 , l ) ˆ 3 . ∗ Dl ( l , : ) ; temp = h l i n i t i a l ( 2 , l ) ; h l i n i t i a l ( 2 , l ) = Beta h ∗ Alpha h ’ / ( Beta h ∗ Beta h ’ ) ; h l i n i t i a l (2 , l ) = conj ( h l i n i t i a l (2 , l ) ) ; h l i n i t i a l ( 1 , l ) = temp ; % update f l A l p h a f = Alpha h − h l i n i t i a l ( 2 , l ) . ∗ A l ( l , : ) ; Beta f = h l i n i t i a l (2 , l ) .∗ Bl ( l , : ) ; E t a f = h l i n i t i a l ( 2 , l ) . ∗ Cl ( l , : ) ; T h e t a f = h l i n i t i a l ( 2 , l ) . ∗ Dl ( l , : ) ; a = Theta f ∗ Theta f ’ ∗ 6; b = ( Eta f ∗ Theta f ’ + Theta f ∗ Eta f ’ ) ∗ 5; c = ( Beta f ∗ Theta f ’ + Eta f ∗ Eta f ’ + Theta f ∗ Beta f ’ ) ∗ 4; d = ( Beta f ∗ Eta f ’ + Eta f ∗ Beta f ’ − Alpha f ∗ Theta f ’ − Theta f ∗ Alpha f ’ ) ∗ 3; e = ( Beta f ∗ Beta f ’ − Alpha f ∗ Eta f ’ − Eta f ∗ Alpha f ’ ) ∗ 2; f = −( A l p h a f ∗ B e t a f ’ + B e t a f ∗ A l p h a f ’ ) ; ppp = [ a b c d e f ] ; r r r = r o o t s ( ppp ) ; temp = f l i n i t i a l ( 2 , l ) ; f l i n i t i a l (2 , l ) = r r r ( length ( r r r ) ) ; f l i n i t i a l ( 1 , l ) = temp ; end i f ( max ( abs ( h l i n i t i a l ( 2 , : ) − h l i n i t i a l ( 1 , : ) )  . / abs ( h l i n i t i a l ( 1 , : ) ) ) < d e l t a h ) . . .  && ( max ( abs ( f l i n i t i a l ( 2 , : ) − f l i n i t i a l ( 1 , : ) )  . / abs ( f l i n i t i a l ( 1 , : ) ) ) < d e l t a f )  count = i t e r a t i o n ; i t e r a t i o n = 0; else i t e r a t i o n = i t e r a t i o n − 1; end end hl hat = h l i n i t i a l (2 , : ) ; f l h a t = f l i n i t i a l (2 , : ) ; end  85  Appendix E. Matlab Code f u n c t i o n [ h l h a t , f l h a t , count ] = c h a n n e l e s t A 2 ( N f f t , L , p , Al , Bl , Cl , Dl , Y , d e l t a h , d e l t a f , h l i n i t i a l ,  ...  fl initial )  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%  %%  %% F i l e : c h a n n e l e s t A 2 .m  %%  %%  %%  %% F u n c t i o n : E s t i m a t e t h e channel parameters ( h l & f l ) u s i n g A l g o r i t h m 2.%% %%  %%  %% Parameters :  %%  %%  N f f t => Size o f FFT  %%  %%  L => L−path channel  %%  %%  p => t h e o r d e r o f T a y l o r S e r i e s Expansion  %%  %%  Al , Bl , Cl , Dl => t h e c o e f f i c i e n t s o f t h e a p p r o x i m a t i o n  %%  %%  Y => t h e r e c e i v e d s i g n a l  %%  %%  d e l t a h => t h e t h r e s h o l d o f h l ( s t o p s t h e i t e r a t i o n )  %%  %%  d e l t a f => t h e t h r e s h o l d o f f l  %%  %%  h l i n i t i a l => t h e i n i t i a l v a l u e o f e s t i m a t e ’ h l h a t ’  %%  %%  f l i n i t i a l => t h e i n i t i a l v a l u e o f e s t i m a t e ’ f l h a t ’  %%  %%  %%  %% Outputs : The esimates ’ h l h a t ’ ,  ’ f l h a t ’ and t h e i t e r a t i o n  ’ count ’  %%  %% %%  %% N o t i c e : The o r d e r o f accuracy depends on N f f t and p . We have  %%  %%  d i f f e r e n t accuracy f o r d i f f e r e n t path . The i t e r a t i o n t i m e  %%  %%  depends on t h e i n i t i a l v a l u e and t h e t h r e s h o l d .  %%  %%  %%  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  %%% I n i t i a l i z e h l and f l %%% i f n a r g i n == 10 YYY = Y . ’ ; AAA = A l . ’ ; BBB = B l . ’ ; temp H = i n v ( A l ∗ AAA) ∗ A l ∗ YYY ; temp F = i n v ( B l ∗ BBB) ∗ B l ∗ (YYY − AAA ∗ temp H ) . / temp H ; h l i n i t i a l = [ temp H temp H ] ; f l i n i t i a l = [ temp F temp F ] ; hl initial = hl initial . ’ fl initial = fl initial . ’ end eta = 4 / 3; count = 0 ; h l h a t = zeros ( 1 , L ) ;  86  Appendix E. Matlab Code f l h a t = zeros ( 1 , L ) ;  switch p case 2  % second o r d e r a p p r o x i m a t i o n  i t e r a t i o n = 1000; while ( i t e r a t i o n ) for l = 1 : L Alpha h = zeros ( 1 , N f f t ) ; A l p h a f = zeros ( 1 , N f f t ) ; if  l == 1 for k = 2 : L Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k ,  :)...  + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) ) ; end e l s e i f l == L for k = 1 : L − 1 Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k ,  :)...  + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) ) ; end else for k = 1 : l − 1 Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k ,  :)...  + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) ) ; end for k = l + 1 : L Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k ,  :)...  + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) ) ; end end % update h l Alpha h = Y − Alpha h ; Beta h = A l ( l , : ) + f l i n i t i a l ( 2 , l ) . ∗ B l ( l , : ) + f l i n i t i a l ( 2 , l ) ˆ 2 . ∗ Cl ( l , : ) ;  temp1 = Beta h ∗ Alpha h ’ / ( Beta h ∗ Beta h ’ ) ;  %  h ( i +1)  temp1 = c o n j ( temp1 ) ; temp2 = e t a ∗ temp1 + ( 1 − e t a ) ∗ h l i n i t i a l ( 2 , l ) ; temp = h l i n i t i a l ( 2 , l ) ; h l i n i t i a l ( 1 , l ) = temp ; h l i n i t i a l ( 2 , l ) = temp2 ; % update f l A l p h a f = Alpha h − h l i n i t i a l ( 2 , l ) . ∗ A l ( l , : ) ; Beta f = h l i n i t i a l (2 , l ) .∗ Bl ( l , : ) ; E t a f = h l i n i t i a l ( 2 , l ) . ∗ Cl ( l , : ) ;  87  Appendix E. Matlab Code a = Eta f ∗ Eta f ’ ∗ 4; b = ( Beta f ∗ Eta f ’ + Eta f ∗ Beta f ’ ) ∗ 3; c = ( Beta f ∗ Beta f ’ − Alpha f ∗ Eta f ’ − Eta f ∗ Alpha f ’ ) ∗ 2; d = −( A l p h a f ∗ B e t a f ’ + B e t a f ∗ A l p h a f ’ ) ; ppp = [ a b c d ] ; r r r = r o o t s ( ppp ) ; temp1 = r r r ( l e n g t h ( r r r ) ) ;  %  f ( i +1)  temp2 = e t a ∗ temp1 + ( 1 − e t a ) ∗ h l i n i t i a l ( 2 , l ) ; temp = f l i n i t i a l ( 2 , l ) ; f l i n i t i a l ( 1 , l ) = temp ; f l i n i t i a l ( 2 , l ) = temp1 ; end i f ( max ( abs ( h l i n i t i a l ( 2 , : ) − h l i n i t i a l ( 1 , : ) )  . / abs ( h l i n i t i a l ( 1 , : ) ) ) < d e l t a h ) . . .  && ( max ( abs ( f l i n i t i a l ( 2 , : ) − f l i n i t i a l ( 1 , : ) )  . / abs ( f l i n i t i a l ( 1 , : ) ) ) < d e l t a f )  count = i t e r a t i o n ; i t e r a t i o n = 0; else i t e r a t i o n = i t e r a t i o n − 1; end end hl hat = h l i n i t i a l (2 , : ) ; f l h a t = f l i n i t i a l (2 , : ) ;  case 3  % t h i r d order approximation  i t e r a t i o n = 1000; while ( i t e r a t i o n ) for l = 1 : L Alpha h = zeros ( 1 , N f f t ) ; A l p h a f = zeros ( 1 , N f f t ) ; if  l == 1 for k = 2 : L Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k , : ) + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) . . . + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 3 . ∗ Dl ( k , : ) ) ; end  e l s e i f l == L for k = 1 : L − 1 Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k , : ) + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) . . . + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 3 . ∗ Dl ( k , : ) ) ; end else for k = 1 : l − 1 Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k , : ) + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) . . . + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 3 . ∗ Dl ( k , : ) ) ;  88  Appendix E. Matlab Code end for k = l + 1 : L Alpha h = Alpha h + h l i n i t i a l ( 2 , k ) . ∗ ( A l ( k , : ) + f l i n i t i a l ( 2 , k ) . ∗ B l ( k , : ) . . . + f l i n i t i a l ( 2 , k ) ˆ 2 . ∗ Cl ( k , : ) + f l i n i t i a l ( 2 , k ) ˆ 3 . ∗ Dl ( k , : ) ) ; end end % update h l Alpha h = Y − Alpha h ; Beta h = A l ( l , : ) + f l i n i t i a l ( 2 , l ) . ∗ B l ( l , : ) + f l i n i t i a l ( 2 , l ) ˆ 2 . ∗ Cl ( l ,  :)...  + f l i n i t i a l ( 2 , l ) ˆ 3 . ∗ Dl ( l , : ) ;  temp1 = Beta h ∗ Alpha h ’ / ( Beta h ∗ Beta h ’ ) ;  %  h ( i +1)  temp1 = c o n j ( temp1 ) ; temp2 = e t a ∗ temp1 + ( 1 − e t a ) ∗ h l i n i t i a l ( 2 , l ) ; temp = h l i n i t i a l ( 2 , l ) ; h l i n i t i a l ( 1 , l ) = temp ; h l i n i t i a l ( 2 , l ) = temp2 ; % update f l A l p h a f = Alpha h − h l i n i t i a l ( 2 , l ) . ∗ A l ( l , : ) ; Beta f = h l i n i t i a l (2 , l ) .∗ Bl ( l , : ) ; E t a f = h l i n i t i a l ( 2 , l ) . ∗ Cl ( l , : ) ; T h e t a f = h l i n i t i a l ( 2 , l ) . ∗ Dl ( l , : ) ; a = Theta f ∗ Theta f ’ ∗ 6; b = ( Eta f ∗ Theta f ’ + Theta f ∗ Eta f ’ ) ∗ 5; c = ( Beta f ∗ Theta f ’ + Eta f ∗ Eta f ’ + Theta f ∗ Beta f ’ ) ∗ 4; d = ( Beta f ∗ Eta f ’ + Eta f ∗ Beta f ’ − Alpha f ∗ Theta f ’ − Theta f ∗ Alpha f ’ ) ∗ 3; e = ( Beta f ∗ Beta f ’ − Alpha f ∗ Eta f ’ − Eta f ∗ Alpha f ’ ) ∗ 2; f = −( A l p h a f ∗ B e t a f ’ + B e t a f ∗ A l p h a f ’ ) ; ppp = [ a b c d e f ] ; r r r = r o o t s ( ppp ) ; temp1 = r r r ( l e n g t h ( r r r ) ) ;  %  f ( i +1)  temp2 = e t a ∗ temp1 + ( 1 − e t a ) ∗ h l i n i t i a l ( 2 , l ) ; temp = f l i n i t i a l ( 2 , l ) ; f l i n i t i a l ( 1 , l ) = temp ; f l i n i t i a l ( 2 , l ) = temp1 ; end i f ( max ( abs ( h l i n i t i a l ( 2 , : ) − h l i n i t i a l ( 1 , : ) )  . / abs ( h l i n i t i a l ( 1 , : ) ) ) < d e l t a h ) . . .  && ( max ( abs ( f l i n i t i a l ( 2 , : ) − f l i n i t i a l ( 1 , : ) )  . / abs ( f l i n i t i a l ( 1 , : ) ) ) < d e l t a f )  count = i t e r a t i o n ; i t e r a t i o n = 0; else i t e r a t i o n = i t e r a t i o n − 1; end  89  Appendix E. Matlab Code end hl hat = h l i n i t i a l (2 , : ) ; f l h a t = f l i n i t i a l (2 , : ) ; end  90  

Cite

Citation Scheme:

    

Usage Statistics

Country Views Downloads
United States 24 0
China 17 4
India 7 0
Canada 6 3
Finland 6 0
France 6 0
Germany 4 4
Vietnam 3 1
Iran 2 1
Algeria 1 0
Mauritius 1 0
Pakistan 1 0
City Views Downloads
Unknown 18 5
Beijing 17 0
Ashburn 12 0
Washington 7 0
Tampere 4 0
Vancouver 3 2
Hanoi 3 1
Kelowna 3 1
Chennai 3 0
University Park 2 0
Rose Hill 1 0
Mountain View 1 0
Dallas 1 0

{[{ mDataHeader[type] }]} {[{ month[type] }]} {[{ tData[type] }]}
Download Stats

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0067250/manifest

Comment

Related Items