Hybrid Spherical- and Planar-Wave Channel Modeling and Estimation for Terahertz Integrated UM-MIMO and IRS Systems

In the literature, mainly two categories of MIMO channel models are considered, namely, the spherical-wave model (SWM) and the planar-wave model (PWM), which are effective in addressing the near-field and far-field effects, respectively [17, 18]. As an improvement to PWM and SWM, we proposed a hybrid spherical- and planar-wave channel model (HSPM) for THz UM-MIMO systems in [13] , which accounted for PWM within the subarray and SWM among subarrays. Compared to the PWM and SWM, the HSPM is more effective by deploying a few channel parameters to achieve high accuracy in the near-field. In the IRS assisted communication systems, an alternative physically feasible Rayleigh fading model was proposed in [12] under the far-field assumption. By taking both near-field effect and IRS into consideration, the authors in [19] considered the SWM for THz integrated IRS and UM-MIMO systems. However, the SWM suffers from high complexity with the massive number of elements in the UM-MIMO and IRS [13]. To date, an effective model addressing the near-field effect in UM-MIMO and IRS systems is still required.

In the IRS systems, the channel analysis mainly focuses on sum rate, power gain, spectral efficiency (SE), and energy efficiency (EE). In microwave systems, the authors in [20] characterized the capacity limit by jointly optimizing the IRS reflection coefficients and the MIMO transmit covariance matrix. The distribution and the outage probability of the sum rate were derived in [21], by considering the SWM of the LoS and PWM of the NLoS. A closed-form expression of the power gain was derived in [22], and the near-field and far-field behaviors were analyzed. At higher frequencies, the ergodic capacity under the Saleh-Valenzuela model was derived and optimized in [23], while the SE and EE are analyzed in [19]. As a critical metric to assess the spatial-multiplexing capability of the channel, the channel rank analysis has been conducted in the THz UM-MIMO systems. To enhance the limited spatial multiplexing in the THz UM-MIMO systems, a widely-spaced multi-subarray (WSMS) structure with enlarged subarray spacing was proposed in [24], where the channel rank can be improved by a factor equal to the number of subarrays. However, the rank analysis in the THz integrated UM-MIMO and IRS systems are still lacking in the literature.

CE for IRS assisted MIMO systems has been explored in the literature [25, 26, 27, 28, 29, 30, 31, 32, 33, 34], which can be categorized into two main categories, namely, estimation of the segmented channels from user equipment (UE) to IRS and IRS to base station (BS), and estimation of the UE-IRS-BS cascaded channel. On one hand, since the passive IRS lacks signal processing capability, it is hard to directly separate each channel segment. Thus, the segmented CE schemes often require special hardware design, e.g., inserting active IRS elements or using full-duplex equipment, both of which however increase the hardware cost [25, 26, 27, 28]. In [25] and [26], a few IRS elements were activated during the pilot reception. The deep-learning tool was then assisted for CE with considerable estimation accuracy. By deploying a full-duplex operated BS, a two timescale CE method was proposed in [27]. The segmented CE problem was formulated as a matrix factorization problem and solved in [28], which can be operated with purely passive IRS. However, this scheme does not address the near-field effect.

On the other hand, since most precoding designs are based on the knowledge of the cascaded channel, the estimation of which has been explored in most existing schemes [29, 30, 31, 32, 33, 34]. In [29], a two-stage atomic norm minimization problem was formulated, by which the super-resolution channel parameter estimation was conducted to efficiently obtain the channel-state-information. Theoretical analysis of the required pilot overhead and a universal CE framework were proposed in [30], which are effective in guiding the design of training and CE. However, all of them are limited to be applicable with fully digital MIMO structures. By exploiting the channel sparsity in the mmWave and THz bands, compressive sensing (CS) based CE methods in hybrid MIMO systems were explored in [31, 32, 33, 34]. These schemes deploy the spatial discrete Fourier transform (DFT) based on-grid codebook to sparsely represent the channel, which is beneficial in achieving reduced training overhead. On the downside, the near-field effect was not incorporated in the DFT codebook, which results in limited estimation accuracy of these schemes in the near-field region of the THz multi-antenna systems [13].

The PWM suitable for the far-field propagation region can be denoted as [18]

where $\alpha_{p}$ represents the complex gain of the $p^{\rm th}$ propagation path, $p=1,...,N_{p}$, with $N_{p}$ denoting the total number of paths. The array steering vectors at Rx and Tx are denoted as $\mathbf{a}_{rp}=\mathbf{a}_{N_{r}}(\psi_{rpx},\psi_{rpz})\in\mathbb{C}^{N_{r}}$ and $\mathbf{a}_{tp}=\mathbf{a}_{N_{t}}(\psi_{tpx},\psi_{tpz})\in\mathbb{C}^{N_{t}}$, respectively. Without loss of generality, by considering an $N$ element planar-shaped array on the x-z plane with physical angle pair $(\theta,\phi)$, the array steering vector $\mathbf{a}_{N}(\psi_{x},\psi_{z})\in\mathbb{C}^{N}$ can be expressed as

where $\psi_{n}=\frac{d_{n_{x}}}{\lambda}\psi_{x}+\frac{d_{n_{z}}}{\lambda}\psi_{z}$, $\psi_{x}={\rm sin}\theta{\rm cos}\phi$, $\psi_{z}={\rm sin}\phi$ denotes the virtual angles, $d_{n_{x}}$ and $d_{n_{z}}$ stand for the distances between the $n^{\rm th}$ antenna to the first antenna on x- and z-axis, respectively.

The SWM is universally applicable to different communication distances, which individually calculates the channel responses of all antenna pairs between Tx and Rx to obtain the ground-truth channel. Due to the high complexity, the SWM is usually deployed in the near-field region, where the PWM becomes inaccurate. By denoting the communication distance of the $p^{\rm th}$ propagation path from the $n_{t}^{\rm th}$ transmitted antenna to the $n_{r}^{\rm th}$ received antenna as $D^{n_{t}n_{r}}_{p}$, the channel response of each antenna pair can be depicted as [18]

where $\mathbf{H}_{\rm S}\in\mathbb{C}^{N_{r}\times N_{t}}$ denotes the SWM channel matrix, $\alpha^{n_{r}n_{t}}_{p}$ represents the complex path gain.

The HSPM accounts for the PWM within one subarray and the SWM among subarrays, which can be denoted as [13]

where $D^{k_{r}k_{t}}_{p}$ stands for the distance between the $k_{r}^{\rm th}$ received and $k_{t}^{\rm th}$ transmitted subarray. The array steering vectors of the $p^{\rm th}$ path for the corresponding subarray pairs are denoted as $\mathbf{a}^{k_{r}k_{t}}_{rp}=\mathbf{a}_{N_{ar}}(\psi^{k_{r}k_{t}}_{rpx},\psi^{k_{r}k_{t}}_{rpz})$, and $\mathbf{a}_{tp}^{k_{r}k_{t}}=\mathbf{a}_{N_{at}}(\psi^{k_{r}k_{t}}_{tpx},\psi^{k_{r}k_{t}}_{tpz})$, respectively, which have similar forms as (4). The virtual angles $\psi^{k_{r}k_{t}}_{rpx}={\rm sin}\theta^{k_{r}k_{t}}_{rp}{\rm cos}\phi^{k_{r}k_{t}}_{rp}$, $\psi^{k_{r}k_{t}}_{rpz}={\rm sin}\phi^{k_{r}k_{t}}_{rp}$, $\psi^{k_{r}k_{t}}_{tpx}={\rm sin}\theta^{k_{r}k_{t}}_{tp}{\rm cos}\phi^{k_{r}k_{t}}_{tp}$, $\psi^{k_{r}k_{t}}_{tpz}={\rm sin}\phi^{k_{r}k_{t}}_{tp}$, where $(\theta^{k_{r}k_{t}}_{rp},\phi^{k_{r}k_{t}}_{rp})$ and $(\theta^{k_{r}k_{t}}_{tp},\phi^{k_{r}k_{t}}_{tp})$ stand for the azimuth and elevation angle pairs at Rx and Tx, respectively. Moreover, $N_{ar}$ and $N_{at}$ depict the number of antennas on the subarrays at Rx and Tx, respectively. We point out that the PWM and SWM are two special cases of HSPM when $K_{t}=K_{r}=1$ and $K_{t}=N_{t}$, $K_{r}=N_{r}$. In addition, the HSPM is accurate and can be adopted when the communication distance is smaller than the Rayleigh distance.

The rank of $\mathbf{H}_{\rm HSPM}$ in (6) satisfies $\min\{K_{r}N_{p},K_{t}N_{p},N_{r},N_{t}\}\leq{\rm Rank}(\mathbf{H}_{\rm HSPM})\leq\min\{K_{r}K_{t}N_{p},N_{r},N_{t}\}$.

Proof: The dimension of the channel matrix $\mathbf{H}_{\rm HSPM}$ in (6) is $N_{r}\times N_{t}$, the maximum rank of the channel is $\min\{N_{r},N_{t}\}$. To prove the right-hand side inequality, we first consider for a fixed transmit subarray $k_{t}$ and propagation path $p$, elements in the set of array steering vectors $\left\{\mathbf{a}_{tp}^{1k_{t}},\dots,\mathbf{a}_{tp}^{K_{r}k_{t}}\right\}$ are linearly independent. Due to the distinguishability of propagation paths, the angles of different paths are different. Therefore, by fixing transmit and receive subarrays as $k_{r}$ and $k_{t}$, respectively, $\mathbf{a}_{tp}^{k_{r}k_{t}}$ for different propagation paths $p=1,\dots,N_{p}$ are linearly independent.

We consider that the $n_{r}^{\rm th}$ received antenna is the $n_{ar}^{\rm th}$ element on the $k_{r}^{\rm th}$ received subarray. Therefore, the $n_{r}^{\rm th}$ row of the HSPM channel in (6) $\mathbf{H}_{\rm HSPM}(n_{r},:)$ can be expressed as

Each row of $\mathbf{H}_{\rm HSPM}$ is a linear combination of $K_{r}K_{t}N_{p}$ linearly independent vectors as

where $k_{r}=1,\dots,K_{r}$ and $\mathbf{0}$ is an all-zero vector of dimension $1\times N_{at}$.

However, the angles of different paths to different received subarrays might be the same, leading that vectors in (14) can be linearly dependent, which reduces the rank of the HSPM channel. Thus, there is ${\rm Rank}(\mathbf{H}_{\rm HSPM})\leq\min\{K_{r}K_{t}N_{p},N_{r},N_{t}\}$. To prove the left-hand side inequality, we consider an extreme case. For a fixed propagation path, the angles among different subarray pairs between Tx and Rx are same. In this case, the HSPM equals to the channel model in [24], whose rank has been proved to be equal to $\min\{K_{r}N_{p},K_{t}N_{p},N_{r},N_{t}\}$, which lower bounds the rank of the HSPM. Till here, we have completed the proof for Lemma 1. $\hfill\blacksquare$

For matrices $\mathbf{A}\in\mathbb{C}^{M\times N}$, $\mathbf{B}\in\mathbb{C}^{N\times N}$ and $\mathbf{C}\in\mathbb{C}^{N\times K}$, where $\mathbf{B}$ is a diagonal matrix, and ${\rm rank}(\mathbf{A})=R_{a}$, ${\rm rank}(\mathbf{B})=N$, ${\rm rank}(\mathbf{C})=R_{c}$, we have

where the equality holds when $\mathbf{A}$ and $\mathbf{C}$ are full-rank matrices.

Proof: Since $\mathbf{B}$ is full row rank, we have ${\rm rank}(\mathbf{AB})={\rm rank}(\mathbf{A})=R_{a}$. Then, ${\rm rank}(\mathbf{ABC})\leq\min\{{\rm rank}(\mathbf{AB}),{\rm rank}(\mathbf{C})\}=\min\{R_{a},R_{c}\}$. When $\mathbf{A}$ is a full-rank matrix, we have ${\rm rank}(\mathbf{AB})={\rm rank}(\mathbf{A})=N$. When $\mathbf{C}$ is a full-rank matrix, we have ${\rm rank}(\mathbf{ABC})={\rm rank}(\mathbf{AB})=N$. $\hfill\blacksquare$

Next, we analyze the rank of the cascaded channel. We adopt the PWM in the far-field region, while the SWM and HSPM are deployed for the near-field region, respectively.

In this case, both $\mathbf{H}_{\rm IRS-BS}$ and $\mathbf{H}_{\rm UE-IRS}$ in (2) adopt PWM, whose ranks equal to $N_{p}^{\rm IB}$ and $N_{p}^{\rm UI}$, respectively. By denoting $\mathbf{H}^{\rm cas}$ in (2) as $\mathbf{H}^{\rm cas}_{\rm PWM}$, from Lemma 2, we can state that ${\rm rank}(\mathbf{H}^{\rm cas}_{\rm PWM})\leq\min\{N_{p}^{\rm IU},N_{p}^{\rm BI}\}$. Moreover, when $N_{p}^{\rm IU}=N_{p}^{\rm BI}=N_{p}$, ${\rm rank}(\mathbf{H}^{\rm cas}_{\rm PWM})\leq N_{p}$.

We first consider $\mathbf{H}_{\rm IRS-BS}$ satisfies the far-field condition and adopts the PWM, while $\mathbf{H}_{\rm UE-IRS}$ meets the near-field condition, which deploys SWM or HSPM. Since the number of propagation paths in the THz channel is much smaller than the number of elements in the UM-MIMO and IRS, ${\rm rank}(\mathbf{H}_{\rm IRS-BS})=N^{\rm IB}_{p}<{\rm rank}(\mathbf{H}_{\rm UE-IRS})$. From Lemma 2, we can obtain that ${\rm rank}(\mathbf{H}^{\rm cas})\leq N_{p}^{\rm IB}$. A similar deduction can be drawn when $\mathbf{H}_{\rm IRS-BS}$ meets the near-field condition while $\mathbf{H}_{\rm UE-IRS}$ satisfies the far-field condition. Thus, when $N_{p}^{\rm IU}=N_{p}^{\rm BI}=N_{p}$, there is ${\rm rank}(\mathbf{H}^{\rm cas})\leq N_{p}$.

We denote $\mathbf{H}^{\rm cas}$ as $\mathbf{H}^{\rm cas}_{\rm SWM}$ when both $\mathbf{H}_{\rm IRS-BS}$ and $\mathbf{H}_{\rm UE-IRS}$ meet the near-field condition and adopt SWM. In this case, $\mathbf{H}_{\rm IRS-BS}$, $\overline{\mathbf{P}}$ and $\mathbf{H}_{\rm UE-IRS}$ are full-rank matrices. Based on Lemma 2, we know that ${\rm rank}(\mathbf{H}^{\rm cas}_{\rm SWM})=\min\{M,N_{u},N_{b}\}$. When $N_{u}=N_{b}=M=N$, ${\rm rank}(\mathbf{H}^{\rm cas}_{\rm SWM})=N$. Similarly, we denote $\mathbf{H}^{\rm cas}$ as $\mathbf{H}^{\rm cas}_{\rm HSPM}$ when both $\mathbf{H}_{\rm IRS-BS}$ and $\mathbf{H}_{\rm UE-IRS}$ adopt HSPM. By considering that $K_{m}=K_{u}=K_{b}=K$ and $N_{p}^{\rm IU}=N_{p}^{\rm BI}=N_{p}$, we have ${\rm rank}(\mathbf{H}^{\rm cas}_{\rm HSPM})\leq K^{2}N_{p}$.

From the above analysis, we can state that in the THz integrated UM-MIMO and IRS systems, the total rank of the cascaded channel is limited by the segmented channel with a smaller rank. This suggests that the rank of the cascaded channel is increased only when both segmented channels meet the near-field condition, which inspires us to enlarge the array size and obtain a larger near-field region. It is worth noticing that the above discussions are not dependent on the IRS beamforming matrix $\overline{\mathbf{P}}$. Therefore, we further claim that given fixed segmented channels, the channel rank can not be improved by the IRS.

We will show in Sec. V that the capacity of the THz integrated UM-MIMO and IRS system based on HSPM is close to that based on the ground truth SWM, which reveals the accuracy of the HSPM. In addition, the HSPM possesses lower complexity compared to the SWM [13]. Therefore, we directly adopt the HSPM for both segmented channels during the CE process.

In the literature, the spatial DFT-based on-grid codebook is widely deployed [31, 32, 33, 34]. This codebook treats the entire antenna array as a unit, and considers that the virtual spatial angles $\psi_{x}={\rm sin}\theta{\rm cos}\phi$ and $\psi_{z}={\rm sin}\phi$ are taken form a uniform grid composed of $N_{x}$ and $N_{z}$ points, respectively. $\theta$ and $\phi$ denote the azimuth and elevation angles, while $N_{x}$ and $N_{z}$ refer to the number of antennas on x- and z-axis, respectively. In this way, the channel is sparsely represented as

where $\mathbf{A}_{\rm Dr}\in\mathbb{C}^{N_{r}\times N_{r}}$ and $\mathbf{A}_{\rm Dz}\in\mathbb{C}^{N_{t}\times N_{t}}$ refer to the two-dimensional DFT on-grid codebooks at Rx and Tx, respectively, which hold a similar form, and can be represented as $\mathbf{A}_{\rm D}=\Big{[}\mathbf{a}_{N}(-1,-1)\dots\mathbf{a}_{N}\left(\frac{2(n_{x}-1)}{N_{x}}-1,\frac{2(n_{z}-1)}{N_{z}}-1\right)$ $\dots\mathbf{a}_{N_{x}}\left(\frac{2(N_{x}-1)}{N_{x}}-1,\frac{2(N_{z}-1)}{N_{z}}-1\right)\Big{]}$. The sparse on-grid channel with complex gains on the quantized spatial angles is depicted by $\boldsymbol{\Lambda}_{\rm D}\in\mathbb{C}^{N_{r}\times N_{t}}$.

To assess the performance of the DFT codebook, we first evaluate the sparsity of on-grid channel $\boldsymbol{\Lambda}_{\rm D}$ in (16) in different cases, by considering the HSPM. Moreover, since in practice, there does not exist a grid whose amplitude is strictly equal to 0, we consider that the sparsity of the on-grid channel equals the number of grids whose amplitude is greater than a small value, e.g., 0.01. First, as illustrated in Fig. 2(a), the amplitude of $\boldsymbol{\Lambda}_{\rm D}$ is shown by considering a compact array without enlarging the subarray spacing. The on-grid channel is sparse, the number of grids with an amplitude larger than 0.01 is only 397, which is much smaller than the preset total number of grids, i.e., 262144. By contrast, in Fig. 2(b), the amplitude of the on-grid channel in the WSMS is plotted. The on-grid channel contains 2755 grids with amplitude larger than 0.01. Therefore, the DFT codebook lacks sparsity in representing the HSPM.

We observe that the HSPM in (6) views each subarray as a unit, each block of which is the production of the array steering vectors for the subarrays at Rx and Tx, respectively. Inspired by this, we consider a subarray-based on-grid codebook. At Rx, the virtual spatial angles for each subarray are considered to be taken from fixed $N_{ar}=N_{arx}N_{arz}$ grids, where $N_{arx}$ and $N_{arz}$ refer to the number of elements on x- and z-axis of the subarray at Rx, respectively. The corresponding DFT codebook is expressed as $\mathbf{U}_{\rm Dr}=\Big{[}\mathbf{a}_{N_{ar}}(-1,-1),$ $\dots\!,\mathbf{a}_{N_{ar}}\left(\frac{2(n_{arx}-1)}{N_{arx}}-1,\frac{2(n_{arz}\!-\!1)}{N_{arz}}-1\!\right)\!,$ $\dots,\mathbf{a}_{N_{ar}}\left(\frac{2(N_{arx}-1)}{N_{arx}}-1,\frac{2(N_{arz}-1)}{N_{arz}}-1\right)\Big{]},$ $n_{arx}=1,\dots,N_{arx}$, $n_{arz}=1,\dots,N_{arz}$. We define $\overline{\mathbf{A}}_{\rm r}\in\mathbb{C}^{N_{r}\times N_{r}}$ as the subarray-based codebook at Rx, which deploys $K_{r}$ $\mathbf{U}_{\rm Dr}$ on its diagonal as

The on-grid codebook matrix at Tx $\overline{\mathbf{A}}_{\rm t}\in\mathbb{C}^{N_{t}\times N_{t}}$ is constructed similarly. Therefore, the on-grid representation of the HSPM in (6) based on the subarray-based codebook can be denoted as

where $\overline{\boldsymbol{\Lambda}}\in\mathbb{C}^{N_{r}\times N_{t}}$ is a sparse matrix. If all spatial angles were taken from the grids and not equal to each other, $\overline{\boldsymbol{\Lambda}}$ would contain $K_{r}K_{t}N_{p}$ non-zero elements.

The amplitude of the on-grid channel $\overline{\boldsymbol{\Lambda}}$ in (18) using the proposed codebook is plotted in Fig. 2(c), by considering the same channel as in Fig. 2(b). The number of grids with an amplitude larger than 0.01 is 1609, which is 1164 smaller than that by using the DFT codebook in Fig. 2(b). In addition, to reveal the accuracy of the on-grid channel, we calculate the difference between the real channel ${\mathbf{H}}_{\rm HSPM}$ and the reconstructed channels approximated by the on-grid channel and the codebooks in (16) and (18) as ${\frac{\left\|\mathbf{A}_{\rm Dr}\boldsymbol{\Lambda}_{\rm D}\mathbf{A}_{\rm Dt}^{\rm H}-{\mathbf{H}}_{\rm HSPM}\right\|_{2}^{2}}{{\left\|{\mathbf{H}}_{\rm HSPM}\right\|_{2}^{2}}}}$ and ${\frac{\left\|\overline{\mathbf{A}}_{\rm r}\overline{\boldsymbol{\Lambda}}\overline{\mathbf{A}}_{\rm t}^{\rm H}-{\mathbf{H}}_{\rm HSPM}\right\|_{2}^{2}}{{\left\|{\mathbf{H}}_{\rm HSPM}\right\|_{2}^{2}}}}$, respectively. The approximation error based on the proposed codebook is around 4 dB lower than that based on the DFT codebook. To this end, we state that the proposed codebook is more efficient than the traditional DFT codebook, which possesses higher sparsity and lower approximation error.

We consider an uplink pilot training procedure, as illustrated in Fig. 3, which is conducted in three levels, namely, the BS training, UE training and IRS training, respectively. During the training process, the UE transmits the known pilot signals to the BS via the IRS in $T=T_{b}T_{u}T_{i}$ training slots for uplink CE, where $T_{b}$, $T_{i}$ and $T_{u}$ denote the number of training slots for the BS, IRS and UE, respectively. At the $(b,u,i)^{\rm th}$ slot, the UE deploys the training beamformer $\tilde{\mathbf{F}}_{u}\in\mathbb{C}^{N_{u}\times N_{su}}$ and transmits the pilot signal $\mathbf{s}_{b,u,i}\in\mathbb{C}^{N_{su}}$, $b=1,\dots,T_{b}$, $u=1,\dots,T_{u}$, $i=1,\dots,T_{i}$. In the meantime, the training phase shift vector $\mathbf{p}_{i}\in\mathbb{C}^{M}$ and combiner $\mathbf{W}_{b}\in\mathbb{C}^{N_{b}\times N_{sb}}$ are deployed at the IRS and BS, respectively, to obtain the received signal ${\mathbf{y}}_{b,u,i}\in\mathbb{C}^{N_{sb}}$ at the BS, which can be represented as