An Off-Grid Turbo Channel Estimation Algorithm for Millimeter Wave Communications

The bandwidth shortage has motivated the exploration of the millimeter wave (mmWave) frequency spectrum for future communication networks. To compensate for the severe propagation attenuation in the mmWave band, massive antenna arrays can be adopted at both the transmitter and receiver to provide large array gains via directional beamforming. To achieve such array gains, channel estimation (CE) with high resolution and low latency is of great importance for mmWave communications. However, classic super-resolution subspace CE methods such as multiple signal classification (MUSIC) and estimation of signal parameters via rotation invariant technique (ESPRIT) cannot be applied here due to RF chain constraints. In this paper, an enhanced CE algorithm is developed for the off-grid problem when quantizing the angles of mmWave channel in the spatial domain where off-grid problem refers to the scenario that angles do not lie on the quantization grids with high probability, and it results in power leakage and severe reduction of the CE performance. A new model is first proposed to formulate the off-grid problem. The new model divides the continuously-distributed angle into a quantized discrete grid part, referred to as the integral grid angle, and an offset part, termed fractional off-grid angle. Accordingly, an iterative off-grid turbo CE (IOTCE) algorithm is proposed to renew and upgrade the CE between the integral grid part and the fractional off-grid part under the Turbo principle. By fully exploiting the sparse structure of mmWave channels, the integral grid part is estimated by a soft-decoding based compressed sensing (CS) method called improved turbo compressed channel sensing (ITCCS). It iteratively updates the soft information between the linear minimum mean square error (LMMSE) estimator and the sparsity combiner. Monte Carlo simulations are presented to evaluate the performance of the proposed method, and the results show that it enhances the angle detection resolution greatly.


Introduction
Thanks to the large bandwidth available at millimeter wave (mmWave) frequencies, mmWave communication technology has become a promising technology to meet the experientially increasing demands of future wireless networks [1]. The small wavelength at the mmWave band enables the integration of a massive number of antennas at both the transmitter and receiver. The directional beamforming technique can be used to achieve a sufficient link margin to compensate for the severe propagation attenuation in mmWave bands [2,3]. Since the channel state information, especially the angle of departure (AoD) and angle of arrival (AoA), are essential for the beamforming and coherent detection, CE with high resolution and efficiency is one of the key requirements in mmWave systems. Due to the high cost and high power consumption of mmWave radio frequency (RF) devices, only a very limited number of RF chains can be integrated in the mmWave systems, which results in much less RF chains than antennas. Therefore, classic angle estimation methods with super-resolution, such as multiple signal classification (MUSIC) and estimation of signal parameters via rotation invariant technique (ESPRIT)-based algorithms [4,5], cannot be directly applied in mmWave systems because these methods require the channel observations at each receive antenna in the digital domain. However, in the mmWave systems, the received signals at each antenna are firstly combined in the analog domain and thus no direct output information from each antenna is available for these classic angle estimation methods. Therefore, new CE algorithms with high-resolution should be developed and tailored to mmWave communications.
Recently, many research efforts have been devoted to the CE for the mmWave systems [6-10]. As adopted in IEEE 802. 15.3c standard [6], a polling mechanism is employed to select the best beam vector pair from the known codebooks with p beams at the transmitter and q beams at the receiver. This method consumes pq time slots to achieve angle resolution O(1/p) and O(1/q). In [7], an adaptive compressive sensing (ACS) method has been developed, which iteratively bisections the beam space and tries to find the grid that the path falls in via space match filtering. It needs 2L 3 log 2 (p/L) + 2L 3 log 2 (q/L) time slots to achieve the resolution O(1/p) at the transmitter and O(1/q) at the receiver for an mmWave channel with L paths. This large time slot payload is unacceptable in practical mobile communication systems with massive antenna arrays, especially in the outdoor wireless coverage scenarios. So far, most existing mmWave CE methods, such as [8][9][10], are based on the on-grid channel model, i.e., all path angles are assumed to be located at the quantization grids exactly, which is not the case in practical systems because the angles are usually continuously distributed. When path angles do not lie on the quantization grids, which we refer to as the "off-grid" channel model, a power leakage effect happens and results in the degradation of CE performance [11]. Since the product of the number of simultaneously supported multi-users and the number of the multi-streams depends on the RF chains, it is preferable that higher angle resolution is achieved by advanced channel estimation methods to separate densely distributed users such that spatial multiplexing methods can be employed to support more users [12]. In this "off-grid" case, the gap to the theoretic lower bound is pretty large, and more antennas are required to enhance the resolution of angle estimation at the cost of more hardware and antennas [13]. Therefore, a new method is necessary to enhance the angle detection performance. So far, few works concentrate on the off-grid problem in the mmWave systems with RF constraints. In [14], the continuous angle estimation was formulated as a sparse regularization problem by taking into account the noises in both the observations and the dictionary by modeling the grid points as noisy, which is hard to be applied in CE because the spectral norm of the noise is much bigger than that of the dictionary. In [15], a continuous basis-based compressed sensing (CS) method was proposed to solve the reconstruction problem in sensing applications where the atoms of the dictionary are used to describe a continuous field, such as frequency and angle. However, the first order approximation of Taylor expansion for the array steering vector in the multi-antenna systems is hard to formulate.
In this paper, we study the off-grid problem and propose an enhanced method to improve the resolution of angle estimation in mmWave systems with massive antenna arrays and RF chain constraints. Firstly, we formulate the signal model to describe the off-grid problem in the angle quantization of the mmWave channel, which decomposes each continuously-distributed angle into the integral grid part and fractional off-grid part. Accordingly, an iterative off-grid turbo CE (IOTCE) algorithm is developed to renew and upgrade the CE between the integral grid and fractional off-grid angle information under the Turbo principle. By fully exploiting the sparse structure of mmWave channels, a soft-decoding based CS method, termed improved turbo compressed channel sensing (ITCCS), is developed to estimate the integral grid angle information, which iteratively updates the soft information between the linear minimum mean square error (LMMSE) estimator and the sparsity combiner. Simulation results show that the proposed CE method can achieve a higher resolution than the existing ones.
The rest of the paper is organized as follows. The mmWave system model is introduced in Section 2. Then, a new off-grid formulation is specified in Section 3 and the proposed IOTCE algorithm is detailed in Section 4. Section 5 presents the simulation results and Section 6 concludes the paper.
Notations: Fonts A, a, a and A denote a matrix, a column vector, a scalar and a set, respectively. A T and A H are the transpose and Hermitian of A. 1 N ∈ C N×1 is an all-ones vector and I is an identity matrix. e N (i) ∈ C N×1 is a unit vector with 1 at its i-th entry and 0 elsewhere. The operator ⊗ are used to denote the Kronecker product. Vect(A) is the matrix vectorization operator that transforms a matrix into a vector by stacking its column vectors, and Diag(a) stands for a diagonal matrix with the entries of a on its diagonal. [A] R,P consists of row vectors and column vectors of the matrix A whose indices contained in the set R and P, respectively. E(·) represents the expectation. δ(·) is the Dirac delta function and I(·) is the indicator function.

System Model
According to [16], typical mmWave channels feature much stronger line-of-sight (LoS) path than non-LoS (NLoS) paths, and the numerical results in [17] show that approximating the sparse mmWave channel with only one LoS path results in a little capacity loss when employing the steerable directional antennas. Thus, single-path mmWave channel model is used in this paper. When the mmWave channel cannot be approximated as a single-path channel, the proposed algorithm cannot be applied directly. New algorithms need to be developed to deal with independently but not identically distributed multi-paths. Consider a downlink mmWave system consisting of a base station (BS) and a mobile station (MS), where the BS is configured with N B antennas and K B RF chains with K B N B while the MS is configured with N M antennas and K M RF chains with K M N M , as shown in Figure 1. The BS applies the beamforming matrix W B ∈ C N B ×K B and the MS employs the matrix W M ∈ C N M ×K M to combine the received signal y ∈ C K M ×1 . The received signal at the MS can be written as: where H ∈ C N M ×N B is the channel matrix, x ∈ C K B ×1 represents the transmitted pilot signal, and z ∈ C K N ×1 is the Gaussian noise with E(zz H ) = σ 2 I. In this paper, we also assume: • uniform linear arrays (ULAs) with half-wavelength spacing are deployed at both the BS and MS; • far-field scattering and block-fading are held, which means the signal waves arrive at different antennas with the same fading amplitudes but distinct phases and the channel fadings are kept constant during the CE precedure.
Under these assumptions, the single-path channel between the BS and MS can be expressed as: where β is the fading weight, θ B ∈ [− π 2 , π 2 ), and θ M ∈ [− π 2 , π 2 ) are, respectively, AoD and AoA. a B and a M are the array steering vector at the BS and the combining vector at the MS, respectively.
, a B and a M can be written as: The (m, n)-th element of the channel can be rewritten as: By discretizing Equation (5) at the angle period [− 1 2 , 1 2 ), where δ M and δ B are the discretization errors corresponding to the φ M and φ B , respectively. The on-grid algorithms, such as [6-10], can estimate the AoA/AoD effectively, by ignoring the discretization errors and assuming that the AoA/AoD are taken from a uniform grid of N points, where N stands for the number of antennas and φ stands for the angle. In practice, the number of antennas is limited and the AoA/AoD are actually continuous, which causes the leakage of power, and, in turn, leads to the degradation of CE performance.

Off-Grid Channel Formulation
In order to suppress the negative effect caused by the off-grid problem, a specific and operable virtual channel representation model is introduced. From this virtual model, an enhanced off-grid angle estimation model is developed.
The continuously-distributed AoA can be written as: where φ M,k = k − 1 is the nearest k-th (k = 1, 2, ..., N M ) discrete quantization grid of φ M , which is referred to as the integral grid angle, and |δ M | < 1 2 stands for the deviation between φ M and φ M,k , termed fractional off-grid angle. Similarly, the discrete form of the AoD is formulated as: where φ B,i = i − 1 (i = 1, 2, ..., N B ) and |δ B | < 1 2 are, respectively, the integral grid part and fractional off-grid part of the angle. Then, the fractional off-grid angle vectors d M ∈ C N M ×1 and d B ∈ C N B ×1 can be rewritten as: Accordingly, the channel matrix in Equation (2) can be reformulated as: where  [18], H D can be rewritten as: where F B ∈ C N B ×N B is the discrete Fourier transform (DFT) matrix whose (m B , n B )-th entry is Then, the virtual representation model of the mmWave channel with continuously-distributed angle can be written as: where D M and D B are unitary matrices. It is interesting that the channel matrix H can be regarded as a two-dimensional (2D) orthogonal transform from the spatial domain to the beamforming domain.
More specifically, the 2D orthogonal transform becomes the 2D DFT when the AoA/AoD locate at the quantization grids exactly.

The Off-Grid Turbo Channel Estimation Algorithm
In order to facilitate the operation of the proposed CE algorithm, the codebook designed in [17] is employed to divide the AoA/AoD estimation into two steps to estimate AoA and AoD separately at the BS and MS in adjacent time slots. Moreover, it can also be adopted here to estimate the fractional off-grid parts of AoA and AoD. Without loss of generality, the codebook design for AoA estimation at the MS is introduced as an example, and the same process can be employed at the BS to estimate the AoD. The process of training mode design is summarized as follows: In this way, the received signal at the MS can be written as: Stacking y r into a vector, we have MK T ∈ C RK×N , with N and K being the number of antennas and RF chains, respectively. F ∈ C N×N is a DFT matrix, and g ∈ C N×1 is the integral virtual channel vector only containing the AoA of the LoS path seen by the MS. D ∈ C N×N is the fractional off-grid angle matrix corresponding to the AoA. Apparently, W is a permutation matrix created by choosing K rows from an identity matrix I. From Equation (17), the CE is accomplished by estimating the integral virtual channel, which, in turn, defines integral grid angle, and the fractional off-grid angle matrix, which stands for the fractional off-grid angle. As shown in Figure 2, an iterative off-grid turbo CE (IOTCE) algorithm is developed, which consists of two components: the integral grid angle estimator and the fractional off-grid angle estimator.

The Integral Grid Angle Estimator
Equation (17) degenerates into a CS recovery algorithm with W being the sampling matrix and D × F being the transform matrix when the fractional off-grid angle δ is given. As the application of the TCCS algorithm in [19], the ITCCS algorithm is adopted to iteratively estimate the integral virtual channel vector g from the Equation (17). g ∈ C N×1 is a sparse signal to be estimated with L non-zero elements, following the Bernoulli-Gaussian distribution with λ = L N being the sparsity. Let f Fg and q D f , and Equation (17) can be rewritten as: As illustrated in Figure 3, the ITCCS algorithm consists of two modules: module A, the LMMSE estimator, and module B, the sparsity combiner. The LMMSE estimator provides an estimate of q by the linear minimum mean square estimation (LMMSE) method, from which an estimate of g can be obtained through constraints. Using the probability density estimation method, the sparsity combiner refines the estimate by exploiting the Bernoulli-Gaussian distribution. The specific realization of the algorithm can be seen in [19]. Therefore, by fixing δ and using the constraints Equations (9) and (10), we get the corresponding fractional off-grid angle vector d(δ) and matrix D, to estimate the integral virtual channel vector g at the sparsity λ = L N (L = 1) using the ITCCS algorithm, where the value of δ is zero at the first iteration and is fed back from the fractional off-grid angle estimator at the iteration. For the purposes of this article, ITCCS (W T , Diag(d t ), 1 N , σ 2 ) is used to make the description of the output, which is the estimate of g, and the parameters in brackets present the input of the algorithm.

The Fractional Off-Grid Angle Estimator
By fixing the integral virtual channel g, and as a result of the equivalence of an elementary row operation and premultiplication by a permutation matrix according to the theory of linear algebra, the Equation (17) can be rewritten as: where R = {n 1 , n 2 , ..., n K } is a set who is constituted by the row indices that are selected from I to W and the operator • denotes the element-wise product. Equation (19) can be solved by the least square (LS) approach. According to the proposed iteration, the IOTCE algorithm is illustrated in Algorithm 1, which presents a structure of a nested loop. In addition, the corresponding explanation of parameters in the IOTCE algorithm is listed in Table 1. The inner loop iterates the integral grid angle between sparsity and linearity, while the integral part and the fractional part are updated in turn during each iteration of the outer loop by following a fixed relationship. As the estimation of integral grid angle improves, the estimation of the fractional off-grid angle enhances via the outer loop. Table 1. Explanation of parameters in the iterative off-grid turbo channel estimation algorithm.

Parameters Explanation
N the number of antennas at receiver g the integral virtual channel vector W the combining matrix at receiver δ the fractional off-grid angle d the corresponding fractional off-grid angle

Simulation Results
In this section, we evaluate the performance of the proposed algorithm in the mmWave communication systems. The simulation system consists of a BS and an MS, and the same ULA with N B = N M = N half-wavelength spaced antennas and K B = K M = K RF chains is configured at both the BS and MS. Typical mmWave channels featuring much stronger LoS paths than NLoS paths are simulated, and a one-path channel model is used to approximate the typical mmWave channel in the low signal-to-noise ratio (SNR) region. The fading weight of the path follows Rayleigh distribution with variance 1, and the continuously-valued AoA/AoD are uniformly-distributed in the range [− π 2 , π 2 ). The simulation results of the average angle estimation error (AAEE), defined as E[|φ −φ|], withφ being the estimate of the true angle φ, are presented in Figure 4, where the on-grid algorithm in [17] and the lower bound of the AAEE for on-grid algorithms, which is 1 N , are also presented for reference in two cases of N = 1024 and N = 512. From the curves in Figure 4, we can observe that the proposed method outperforms the on-grid method and achieves AAEE of 10 −2 at SNR γ = −4 dB. As expected, the AAEE of the proposed algorithm exceeds the lower bound of the on-grid algorithms as SNR increases, which validates the AAEE performance of our algorithm. For example, the gaps between the simulated AAEE and the lower bound are 0.00139 and 0.0007 when γ = 20 dB for N = 512 and N = 1024, respectively. The gap shrinks from 0.00139 to 0.0007 by the increase in N, because the less sparsity (λ = 1 N ) can achieve a the better recovery performance at the given samplings according to the theory of CS [20].   The average probability of integral grid point estimation error (APIEE) performance, which is defined as E[I(φ k =φ k )], withφ k being the estimate of the real integral grid angle φ k , is also simulated, and the simulation results are illustrated in Figure 5. As clearly shown, the bigger N leads to the better APIEE performance at a given K. For example, at the APIEE of 10 −2 , the algorithm with N = 1024 almost enhances by 4 dB with respect to N = 512.

Conclusions
The shortage of spectrum at present will be alleviated by advances in the mmWave communication systems, and directional precoding/beamforming with large antenna arrays appears to be inevitable to support longer outdoor links and to provide sufficient received signal power than before. In this paper, an enhanced CE algorithm is developed for the off-grid problem in the mmWave systems with massive antenna arrays and RF chain constraints. First, we developed an off-grid formulation to catch the quantization problem. Then, an iterative method, which depends on the developed off-grid formulation is proposed to renew and upgrade the CE between the integral and fraction angle information under the Turbo principle. Simulation results show the efficiency and high resolution of the proposed method. For future work, it would be more practical to consider the multi-paths channel with independently but not identically distributed paths, and develop an efficient method to calculate the Log Likelihood Ratio (LLR) information.