I Introduction
To meet the needs of future wireless communications, especially various new services and scenarios that are continually emerging, the fifthgeneration (5G) mobile systems come into being. One objective of 5G is to achieve reliable communication for the scenarios with high Doppler spread, such as highmobility scenarios [39, 8]. Orthogonal frequency division multiplexing (OFDM) technology has been widely used to combat the intersymbol interference (ISI) in timeinvariant channels [6, 33, 21, 7]. However, when there is a high Doppler spread in timevariant channels, which would lead to severe intercarrier interference (ICI), the performance of OFDM degrades significantly[36].
To cope with this problem, orthogonal time frequency space (OTFS) modulation technology was proposed [12, 25, 11]
and attracted much attention due to significant advantages in timevariant channels. OTFS converts the signal in the timefrequency domain into a delayDoppler domain and performs modulation accordingly. In particular, each information symbol in the delayDoppler domain is expressed by a pair of orthogonal basis functions in the timefrequency domain, i.e., full diversity in the timefrequency domain. Furthermore, transmitted symbols experience roughly constant channels even with high Doppler spread. Another advantage is that the implementation of OTFS can be achieved simply as an overlay of the OFDM systems. Farhang
et al. investigated the OTFS system based on OFDM [9] and designed a modem scheme with low complexity, where the cyclic prefix (CP) is added in front of each OFDM symbol in an OTFS symbol.Channel estimation is crucial in OTFS systems. In [26] and [28], trainingbased channel estimation methods were investigated. In [26], Murali et al. used the pseudorandom (PN) sequence as pilots to acquire the channel state information (CSI). In [28], Raviteja et al. arranged the pilots and data in the same OTFS symbol and estimated the delayDoppler channel by a threshold method. After the acquisition of the CSI, the symbol detection can be performed. In [26]
, the Markov chain Monte Carlo (MCMC) sampling was utilized for the lowcomplexity symbol detection. In
[29], an explicit derivation of the inputoutput relationship of OTFS systems was given, and the effect of fractional Doppler caused by insufficient sampling of Dopplerdimension was analyzed. They also developed a message passing (MP) method to detect OTFS symbols. However, when they used the rectangular pulseshaping waveforms, the CP was not considered, which would cause ISI and might make the receiver more complicated. Different from the symbol spaced sampling framework in [29], Ge et al. designed a fractionally sampling scheme for symbol detection [10], where only one CP is inserted for the whole OTFS symbol. Compared with [9], less CP overhead in [28, 29, 10] can improve spectral efficiency. In [37], Wei et al. reduced the channel spread caused by the fractional Doppler and improved the channel estimation performance through the design of transmitter and receiver windows. In [18], the embedded pilotaided CSI acquisition scheme proposed in [28] and the MP algorithm proposed in [29] were extended to multipleinput multipleoutput (MIMO) OTFS systems.Massive MIMO technology has been one of the key technologies of 5G owing to the immense improvement of spectrum and power efficiencies [24, 41, 30, 42, 14, 43, 27, 19, 44]. Considering high Doppler spread scenarios, the integration of massive MIMO and OTFS can further improve the performance. To utilize such benefits, the base station (BS) requires downlink CSI for precoding. In traditional timedivision duplexing (TDD) systems, uplink training can be used to acquire the downlink CSI via leveraging the channel reciprocity [24, 23]. However, when systems work in frequencydivision duplexing (FDD) mode, user terminals need the pilots transmitted by the BS to acquire the CSI and feed it back. Thus the downlink channel estimation is necessary for massive MIMOOTFS systems. In [32], Shen et al. demonstrated the 3dimensional sparsity of the downlink massive MIMOOTFS channel, based on which a 3Dstructured orthogonal matching pursuit (3DSOMP) method was proposed to acquire the downlink CSI. However, the pilot used for the channel estimation is the random one, which is hard to be implemented and consumes much memory for its randomness in the practical system. In [22], Liu et al. perform the downlink channel estimation with the help of channel parameters obtained from uplink training. But downlink CSI acquisition can not be performed independently and must be after the uplink channel estimation. In [31], Shan et al. designed a lowoverhead pilot pattern for the CSI acquisition and equalized the received signal from each angle through the channel information obtained. It is noticed that the fractional Doppler is not considered in [32, 22, 31], and the system models are accurate only when the Doppler frequency of each path of the channel is mapped to an integer tap. It is impractical since the resolution of the Doppler axis is not sufficient. In [20], Li et al. developed a path division multiple access (PDMA) scheme by which different users used scheduled delayDoppler domain grids to communicate with BS simultaneously without interuser interference. However, they only considered TDD systems, where downlink CSI can be obtained by the reciprocity between the uplink and downlink channel. Thus only the uplink channel estimation was discussed, and the proposed scheme is not suitable for FDD systems.
In this work, we consider the fractional Doppler caused by insufficient sampling of Dopplerdimension in practical systems. Moreover, for massive MIMOOTFS systems, we propose a downlink CSI acquisition scheme, including deterministic pilot design and channel estimation algorithm. The contributions are summarized as follows:

In the singleinput singleoutput (SISO) case, we analyze the inputoutput relationship of the OTFS system, where the OFDM modem is chosen as the timefrequency modem. Specifically, we consider that the Doppler frequency of each path is mapped to an integer and a fractional Dopper tap and reveal the influence of the subpaths contained in each dominant path. Moreover, the characteristics of subpaths are used to simplify such a system model. Then we extend it to the massive MIMOOTFS system and establish the downlink CSI acquisition model.

To reconstruct the channel accurately, we design a ZadoffChu (ZC) sequence based deterministic pilot, which can achieve better sensing performance and memory consumption saving than the random pilot used in [32]. We first analyze the relationship between the pilot matrix and the position of pilots at each beam. Then we give two conditions that the deterministic pilots need to satisfy to ensure the low coherence between the pilot matrix columns. Next, ZC sequences are chosen as pilots, whose performance is also discussed.

Based on the modified model, we propose a modified sensing matrix based channel estimation (MSMCE) algorithm to acquire CSI. Due to the slow change of delayDoppler channels, path delays and Doppler frequencies of all dominant paths are extracted from the previous channel estimation result and used to modify the sensing matrix for more accurate CSI acquisition. Then we can utilize the estimated CSI to update the path delays and Doppler frequencies, which can be used in the next channel estimation.
The rest of the paper is organized as follows. In Section II, we investigate models for SISO OTFS systems and massive MIMOOTFS systems. Section III establish a downlink CSI acquisition model. Based on the model, we design a deterministic pilot matrix and propose a MSMCE algorithm. Section IV presents the performance of the proposed CSI acquisition scheme by simulation results, and the paper is concluded in Section V.
Notations: The superscripts and
denote the conjugate and conjugatedtranspose operations, respectively. The uppercase (lowercase) boldface letters denote matrices (column vectors).
denotes the imaginary unit. and denote the norm of and , and denotes the norm of . denotes the smallest integer that is not less than , while denotes the largest integer that is not greater than . is the Hadamard product operator, and is the Kronecker product operator. denotes mod , and denotes . The notation is used for definitions, and is used to indicate equivalence or approximate equivalence. denotes the Dirac delta function. denotes the th row of , wihle denotes the th colomn of . denotes the th element of . denotes the th element of .Ii System Model
In this section, some basic definitions and concepts of OTFS are reviewed first. Then we analyze the inputoutput relationship of the SISO OTFS system, where the OFDM modem is selected as timefrequency modem, and the fractional Doppler is considered. Then we extend it to the case of massive MIMO systems.
Iia SISO OTFS System
For the sake of the subsequent derivation, we first define a lattice in the timefrequency domain as , where (seconds) and (Hz) are sampling intervals of the timedimension and the frequencydimension, respectively, , , . Similarly, the lattice in the delayDoppler domain is defined as where and are sampling intervals of the delaydimension and Dopplerdimension, respectively, and , .
We then arrange a set of quadrature amplitude modulated (QAM) symbols in the delayDoppler domain on the lattice
. The inverse symplectic finite Fourier transform (ISFFT) is utilized at OTFS transmitter to convert
to the symbols in the timefrequency domain as [2](1) 
Next, an OFDM modulator can be used as a timefrequency modulator to convert to a transmitted signal with a transmitted waveform as
(2) 
where is the length of CP, is time duration of an OFDM symbol without CP, and is defined as
(3) 
Note that the CP is added in this step.
After the transmitted signal passing through the multipath timevariant channel, the received signal can be obtained as
(4) 
where is the delay, is the Doppler frequency, is impulse response in the delayDoppler domain[17], and is the additive Gaussian noise. Usually, there are only a few scatterers in the transmission environment. Thus the channel can be represented in a sparse way [29]. We assume that the number of the dominant paths between the transmitter and the receiver is , and each consists of subpaths. Hence, is given by
(5) 
where and are the complex path gain and the Doppler frequency of the th subpath of the th dominant path, respectively. All of subpaths in the th dominant path have the same delay [15]. We define the delay and Doppler taps for subpath as
(6) 
where and are integers and represent the indexes of delay and Doppler taps. The real number , whose value range is , is defined as the fractional Doppler. The fractional delay is not considered since the resolution of the delay axis is sufficient so that each path delay can be mapped to an integer delay tap in typical wideband systems [35].
At the receiver, the crossambiguity function between a received waveform and the received signal is computed by the matched filter as
(7) 
where is defined as
(8) 
Note that the CP is removed in this step. By sampling the crossambiguity function, the received data in the timefrequency domain is given by
(9) 
It is worth noting that when the time duration of CP, i.e., is beyond the maximum path delay of all dominant paths, there is no ISI between OFDM symbols within an OTFS symbol at the receiver, which is different from [28, 29, 10].
Next, the symplectic finite Fourier transform (SFFT) can be utilized to map to the symbols in the delayDoppler domain as
(10) 
We define the complex gain of the timevariant channel on the delay tap at time as
(11) 
where is the system sampling interval, and we define as the complex gain of the th subpath of . Through the above discussion, we give the inputoutput relationship of the SISO OTFS system, as shown in Proposition 1.
Proposition 1: For a SISO OTFS system, when the OFDM modem is used as the timefrequency modem, the inputoutput relationship is expressed as
(12) 
where
(13) 
and is the additive noise in the delayDoppler domain.
Proof: See Appendix.
From (IIA), we have the following two findings. First, is related to , which is the received data position along the delaydimension. Second, for the th dominant path, is affected by all subpaths due to the fractional Doppler. These two points make the system model very complicated, which is not conducive to subsequent analysis. Hence, we utilize the characteristics of subpaths to simplify it. Since all subpaths of one dominant path usually originate from the same scattering cluster, the directions of arrival deviation of these subpaths at the user terminal are usually slight [16]. Therefore, we approximate the Doppler frequency of all subpaths of the th dominant path to the same value , i.e., , . Such an approximation simplifies the model in (IIA) to the follows
(14) 
where is the delayDoppler domain channel, which is defined as
(15) 
Equality in (IIA) holds precisely if there is only one subpath in a dominant path (i.e., ).
It is noticed that when , (IIA) is transformed into
(16) 
which is the same as that in [32] and [22]. Therefore, the model in [32] and [22] can be seen as a special case of our proposed inputoutput relationship when . Moreover, comparing (IIA) with (IIA), we can find that the phase compensation in (IIA) uses the information of each dominant path, which makes the model more accurate in the practical system.
Note that we choose the OFDMbased OTFS system similar to [9] in this paper, which can be achieved simply as an overlay of the widely used OFDM system, instead of the spectral efficient OTFS model in [28, 29, 10]. If fewer CPs are used to improve spectral efficiency as in [28, 29, 10], the inputoutput relationship will have an additional case owing to the ISI between OFDM symbols within an OTFS symbol, and the only difference between the two cases is the exponential term, which implies that the basic idea of using the channel path information to modify the phase compensation matrix of the sensing matrix to improve the channel estimation performance, as shown later, is still valid in the spectral efficient OTFS model.
IiB Massive MIMOOTFS System
We consider the downlink transmission in a massive MIMOOTFS system. The BS deploys antennas and serves singleantenna user terminals. We consider the centralized massive MIMO, where the uniform linear array (ULA) is equipped at the BS, and the antenna spacing is set to half wavelength. Without loss of generality, we focus on one user terminal and disregard the dependency of the channel on user index for simplicity.
Similar to (IIA), the symbols received at the user terminal in the delayDoppler domain can be expressed as
(17) 
where is transmitted symbols in the delayDopplerspace domain, is the delayDopplerspace domain channel, which is defined as [13]
(18) 
where is the angle of departure (AoD) of the th subpath. To exploit the sparsity of the beam domain, the delayDopplerbeam domain channel is obtained by applying normalized discrete Fourier transform (DFT) for the delayDopplerspace domain channel along the spacedimension as [32, 34]
(19) 
where is the beam index, and . From (IIB), we can find that the dominant elements of are distributed in the positions where , and , which means that the channel has the 3D sparsity over the delayDopplerbeam domain [32, 1]. By combining (IIB) and (IIB), the received symbols can be rewritten as
(20) 
where
(21) 
Equality in (IIB) holds precisely if there is only one subpath in a dominant path. Similarly, when , the received symbols are given by
(22) 
According to downlink massive MIMOOTFS system models and the sparsity of the delayDopplerbeam domain channel , the downlink CSI acquisition is actually a sparse signal reconstruction problem [32]. A variety of compressive sensing (CS) algorithms can be used to solve this problem. Next, we will give the model of downlink CSI acquisition and propose the corresponding deterministic pilot design and MSMCE algorithm.
Iii Downlink CSI acquisition for Massive MIMOOTFS Systems
In this section, we first model the downlink CSI acquisition as a sparse signal reconstruction problem. To improve sensing performance, the deterministic pilots are designed, and the pilot overhead is also discussed. Then, we propose a MSMCE algorithm to reconstruct the sparse delayDopplerbeam domain channel and analyze its performance.
Iiia Downlink CSI Acquisition Model
Fig. 1 shows the position of pilots in the delayDoppler domain at one antenna. We consider that the position of pilots in the delayDoppler domain at each transmit antenna is the same and given by
(23) 
where and are the initial position of the pilots in the Doppler domain and delay domain, respectively, and and are the lengths of pilots. The guard intervals (i.e., zero symbols) are demanded to eliminate the interference between the data and pilots.
DelayDopplerbeam domain channels have finite support and along the delaydimension and the Dopplerdimension, respectively [12]. Hence for the delaydimension, the guard intervals should be placed at the beginning and the end of the pilots with the length . However, for the Dopplerdimension, is typically noninteger in practical systems, which means that according to the characteristic of in (IIB), the magnitude at each is not zero and decreases as moves away from . Thus we only consider for and replace the values outside the range with zero, and the guard intervals should be placed at the beginning and the end of the pilots with the length . Therefore, the range of and of the delayDopplerbeam domain channel are limited to and , respectively. Moreover, the data is placed in the position except for the pilot and the guard interval.
The position of received symbols for channel estimation is the same as (23). Therefore, (IIB) can be rewritten as
(24) 
where is the pilot matrix with element of index , where , , , and . In (24), and are the vector forms of received symbols and the additive noise with element and of index . is the vector form of the delayDopplerbeam domain channel with element of index , is the phase compensation matrix with element of index , which is defined as
(25) 
where Del is the tap index set of the delay of all dominant paths, i.e., . is selected from the set Dop which consists of the integer and fractional tap indexes of the Doppler frequency of all dominant paths, i.e., . When , according to (IIB), in (24) is converted to with element of index . And in this case, the model is the same as that in [32].
We can rewrite (24) as
(26) 
where . Thus the downlink CSI acquisition of massive MIMOOTFS systems is converted to a sparse signal reconstruction problem, where is the sparse vector to be recovered, and is the sensing matrix. When , is converted to . Hence we have two kinds of sensing matrix that can be used in CS, and the only difference between them is the phase compensation matrix. The phase compensation matrix is more accurate but requires both path delay and Doppler frequency of each dominant path, which can not be directly obtained. In contrast, the phase compensation matrix is irrelevant to the channel, although it is inaccurate. To combine the advantages of both, the MSMCE algorithm is proposed and will be discussed in detail later. Prior to that, we give the deterministic pilot design in the next subsection.
IiiB Deterministic Pilot Design
In order to recover the sparse vector reliably, the sensing matrix must be designed carefully and satisfies the restricted isometry property (RIP) proposed in [3]
. The RIP implies that the coherence (i.e., inner product) between columns in a sensing matrix should be as small as possible to obtain a good sensing performance. Although the random matrix is usually used as the sensing matrix due to its nearoptimality
[4], it is hard to be implemented and costs a lot of memory consumption for its randomness. Therefore, we focus on the deterministic pilot design for practical communication systems.Since the phase compensation matrix of the sensing matrix is related to the channel, which varies with time, the pilot matrix is what we analyzed exactly. The design of pilot matrix is equivalent to the design of transmitted pilots in the delayDopplerbeam domain with their position in the delayDoppler domain at each beam. Note that the pilots in the delayDopplerspace domain (i.e., the pilots in the delayDoppler domain at each BS transmit antenna) can be obtained using (21), and such a transformation will not affect the position of pilots in the delayDoppler domain.
We first analyze the relationship between the pilot matrix and the pilot position at one beam. As shown in Fig. 2, each dashed box contains pilot symbols, which are used to construct the corresponding column of the pilot matrix in Fig. 2. The pilots at each beam can construct columns totally.
However, we can find that if we directly use the position of pilots shown in Fig. 1, there will be some zero symbols in the pilot matrix, which will affect the coherence between columns of the pilot matrix and make the pilot design difficult. Therefore, based on Fig. 1, we propose a new pilot position design shown in Fig. 3. Specifically, along the Dopplerdimension, the first pilots are added to the end (marked in yellow grids in Fig. 3), and the last pilots are added to the beginning (marked in pink grids in Fig. 3) while the guard intervals are placed at the beginning and the end of pilots with the length . Along the delaydimension, the last pilots are added to the beginning (marked in thickblack grids in Fig. 3) while the guard intervals are placed at the end of the pilots with the length of . The data is placed in the position except for the pilot and the guard interval.
Next, we discuss the pilot design at a given beam and between different beams. For a given beam, we assume that the th column of the pilot matrix (i.e., the black dashed box in Fig. 2) constructed by the pilots at this beam shown in Fig. 3 is defined as
(27) 
where denotes a sequence of length and is cyclically shifted by symbols. According to the designed pilot position and its relationship with the pilot matrix, different columns (i.e., dashed boxes with different colors in Fig. 2) can be expressed as
(28) 
Therefore, according to the property of Kronecker product , the pilot sequence should be orthogonal to its cyclically shifted version to ensure the orthogonality between columns of the pilot matrix at the given beam, which is defined as the orthogonality condition.
The above pilot design can be used at up to beams, where and . The difference between these beams is that should be replaced with , or with , where and . Note that at least one of the two pilot sequences needs to be replaced to ensure the orthogonality between columns. However, for the beam other than these beams, a new pair of pilot sequences should be used. And the coherence of the columns constructed by this new pair of pilot sequences and the columns constructed by other pairs of pilot sequences should be as low as possible to ensure the sensing performance of the entire pilot matrix, which is defined as the low coherence condition.
To sum up, the designed deterministic pilots should satisfy both the orthogonality condition and the low coherence condition. ZC sequence is quite suitable for the proposed deterministic pilot design since it is a cyclic orthogonal sequence with good autocorrelation and low crosscorrelation. We define the ZC sequence cyclically shifted by symbols with root and length as , and its th element can be expressed as
(29) 
where . Note that in the current system should be or . The cyclic orthogonality of the ZC sequence makes it satisfy the orthogonality condition. Meanwhile, when the length of ZC sequence is prime, the magnitude of crosscorrelation between two normalized ZC sequences with different roots is [5]. It means that when we need a new pair of pilot sequences, a pair of ZC sequences with a new pair of roots can be chosen, and they satisfy the low coherence condition. Therefore, based on the ZC sequence, the proposed deterministic pilot design is summarized in Algorithm 1, where the pilots at all beams are rearranged as a matrix of size .
Comments
There are no comments yet.