Cooperative Secrecy Beamforming Scheme Resistant to Multi-eavesdroppers in Dense Heterogeneous Networks

1.   引言

太赫兹(TeraHertz, THz)频谱的大带宽和高数据速率，有望满足无线通信系统中爆炸式增长的数据流量需求；然而，THz的极高传播衰减和分子吸收通常会限制信号的传输距离和覆盖范围^[1]。幸运的是，可重构智能表面(Reconfigurable Intelligent Surface, RIS)能以较低的硬件成本和能耗提高无线通信系统的容量和覆盖率^[2,3]。然而，RIS提供的性能增益很大程度上依赖于信道状态信息(Channel State Information, CSI)的准确性，且RIS的无源性和高维级联信道使得CSI的准确获取具有挑战性^[3]。

针对RIS辅助单用户通信的信道估计方案已取得很多有价值的研究成果^[4-6]。然而，文献[4-6]这类针对RIS辅助单用户通信系统的信道估计方案，直接推广到RIS辅助的多用户通信系统中是低效的^[7]。在已有的RIS辅助多用户通信的宽带信道估计方案中，文献[8]利用基站(Base Station, BS)到RIS的信道高维准静态特性和RIS到用户设备(User Equipment, UE)的信道移动低维特性，提出一种低导频开销的双时间尺度信道估计方案，利用坐标下降算法估计BS-RIS信道，利用最小二乘(Least Square, LS)算法估计RIS-UE和BS-UE信道。文献[9]研究载波频偏(Carrier Frequency Offset, CFO)对RIS辅助的正交频分复用(Orthogonal Frequency Division Multiplexing, OFDM)无线通信系统中信道估计的影响，并提出一种联合CFO和信道冲激响应的信道估计方法。文献[10]提出一种基于深度去噪神经网络的压缩信道估计方法，以减少训练开销。文献[11]将训练信号建模为张量模型，并利用代数算法解决张量补全问题，从而恢复出信道的多径参数。

上述宽带信道估计方案存在以下局限性：首先，文献[8-11]的宽带信道估计方案考虑了宽带信道的频率选择性，但忽略了波束斜视的影响^[12]。在OFDM宽带系统中，波束斜视使阵列响应变得频率相关，导致传统的信道估计方案失效^[13]。其次，文献[8-11]的宽带信道估计方案未利用信道间的强相关性，虽然文献[14]利用了信道间的强相关性，但所提方案是针对窄带系统的。最后，针对RIS辅助的宽带通信系统且考虑波束斜视的研究仅有文献[15]，但文献[15]是针对单用户系统的。因此，针对波束斜视影响下RIS辅助多用户通信的宽带信道估计方案的研究具有非常重要的理论意义与现实价值。

本文研究了波束斜视影响下RIS辅助多用户THz通信的上行信道估计问题，并提出一种低复杂度的稳健两阶段信道估计方案。在第1阶段，利用THz信道稀疏性，将信道估计问题转化为一个非凸目标优化问题，提出一种基于梯度下降法的迭代加权的对数和约束信道估计算法，本算法以对数和函数为约束条件，用优化最小化(Majorization-Minimization, MM)方法得到一个代理目标函数，然后通过梯度下降法优化目标函数，从而准确估计出信道增益、角度信息和时延参数。在第2阶段，利用所有用户在各自的级联信道(UE-RIS-BS)中共享公共信道RIS-BS信道这一特性，将其他用户的级联信道表示为典型用户CSI的缩放版本，利用LS算法进行信道估计。仿真结果表明，本文算法具有较好的性能。

符号说明：$\otimes$为克罗内克积，$\odot$为哈达玛积，$\text{diag}(\cdot)$为对角矩阵，$\text{blkdiag}(\cdot)$为块对角矩阵，${\Vert \cdot\Vert }_{\text{F}}$为弗罗贝尼乌斯范数，$\text{vec}(\cdot)$为列序向量化，$\left\lceil a \right\rceil$为大于或等于$a$的最小整数，$\text{rank}(\cdot)$为矩阵的秩，${(\cdot)}^{†}$为伪逆。

2.   系统模型

本文研究RIS辅助多用户大规模多输入单输出THz无线通信系统中的信道估计问题。在本系统中，BS在RIS的辅助下同时为$U$个单天线UE提供服务。RIS为配备$M = {M_{\text{r}}} \times {M_{\text{o}}}$个无源反射元件的均匀平面阵列，BS为配备${N_{{\text{RF}}}}$个射频链(Radio Frequency, RF)和$N$根天线的均匀线性阵列，RIS上反射元件之间的间距或天线阵列之间的间距为半波长，即$d = \lambda /2$，$\lambda$为载波波长。为对抗多径效应，采用具有$Q$个子载波的OFDM技术。BS和RIS安装位置固定，且考虑低移动性UE。考虑准静态块衰落信道，每个信道在${T_{{\text{all}}}}$个时隙的信道相干块中保持近似恒定^[16]。假设每个相干块中${T_{{\text{up}}}}$个时隙用于上行信道估计，每个时隙包含$P$个导频符号。假设UE-BS之间的直接路径被大的障碍物阻塞^[6]。

在上行信道估计期间，在第$q$个子载波的第$b$个时隙BS的接收信号为

其中，${{\boldsymbol{H}}_{{\text{RU}},u,q}}$为第$u$个UE的第$q$个子载波的UE到RIS的信道矩阵，${{\boldsymbol{H}}_{{\text{BR,}}q}}$为RIS到BS的信道矩阵，${\boldsymbol{\varPhi}} = {\left[ {{\eta _1}{{\text{e}}^{{\text{j}}{\varphi _1}}},{\eta _2}{{\text{e}}^{{\text{j}}{\varphi _2}}}, \cdots ,{\eta _M}{{\text{e}}^{{\text{j}}{\varphi_M}}}} \right]^{\rm T}}$为RIS的相移系数，${\eta _m} \in \left[ {0,1} \right]$和${\varphi _m} \in [0,2\pi )$为第$m$个反射元件的幅度和相移,${{\boldsymbol{w}}_b} \in {\mathbb{C}^{N \times {N_{{\text{RF}}}}}}$为BS第$b$个时隙处的组合矩阵；${{\boldsymbol{s}}_{u,b,q}} = \left[ {{s_{u,b,q,1}},{s_{u,b,q,2}}, \cdots ,{s_{u,b,q,P}}} \right]$为第$u$个UE的第$q$个子载波在第$b$个时隙处发送的导频序列,${{\boldsymbol{n}}_{b,q}}$为加性高斯白噪声矩阵。假设采用正交导频传输策略，即${{\boldsymbol{s}}_{u,b,q}}{\boldsymbol{s}}_{u,b,q}^{\text{H}} = RP$, ${{\boldsymbol{s}}_{u1,b,q}}{\boldsymbol{s}}_{u2,b,q}^{\rm H} = 0,\forall u1 \ne u2$, $R$为发射功率^[15]。因此，在第$q$个子载波的时隙$b$处BS接收到的信号为

其中，${\widetilde {\boldsymbol{y}}_{u,b,q}} \triangleq \dfrac{1}{{RP}}{{\boldsymbol{y}}_{u,b,q}}{\boldsymbol{s}}_{u,b,q}^{\text{H}}$, ${\widetilde {\boldsymbol{n}}_{u,b,q}} \triangleq \dfrac{1}{{RP}}{{\boldsymbol{n}}_{u,b,q}} \cdot {\boldsymbol{s}}_{u,b,q}^{\text{H}}$。

3.   宽带信道模型

假设该系统的传输带宽为$B$，载波频率为${f_{\text{c}}}$，第$q$个子载波的频率为${f_q} = {f_{\text{c}}} + B\left( {q - 1 - \left( {Q - 1} \right)/2} \right)/Q, \forall q \in Q$。本文只考虑BS和RIS之间的视距传输^[3,17]。假设UE和RIS之间存在$L$条路径，$\tau _{u,l,m}^{{\text{RU}}}$表示第$u$个UE到RIS第$m$个反射元件在第$l$条路径上的路径时延，$\tau _{u,m,n}^{{\text{BR}}}$表示RIS的第$m$个反射元件到BS的第$n$根天线的路径时延。因此，第$u$个UE到RIS的第$m$个反射元件的上行链路信道的时域响应为^[13]

其中，${\alpha _{{\text{RU,}}l}}$为UE-RIS信道的第$l$条路径的信道增益。

对式(3)进行傅里叶变换，则第$u$个UE的第$q$个子载波在UE-RIS信道整体频域信道响应为

其中，${\overline \alpha _{_{{\text{RU,}}l}}} \triangleq {\alpha _{{\text{RU,}}l}}{{\text{e}}^{ - {\text{j}}2\pi {f_{\text{c}}}\tau _{u,l}^{{\text{RU}}}}}$为第$l$条路径等效路径增益，且$\tau _{u,l}^{{\text{RU}}} \triangleq \tau _{u,l,0}^{{\text{RU}}}$，${\theta _l}$和${\phi _l}$分别为RIS第$l$路径的到达角(Angle-Of-Arrival, AOA)的仰角和方位角，${{\boldsymbol{a}}_{{\text{RIS}}}}\left( {{\theta _l},{\phi _l},{f_q}} \right)$为RIS的阵列导向矢量矩阵。${{\boldsymbol{a}}_{{\text{RIS}}}}\left( {\theta ,\phi ,{f_q}} \right)$表示为^[13]

同理，第$u$个UE的第$q$个子载波在RIS-BS信道整体频域信道响应为

其中，${\overline \beta _{{\text{BR}}}} \triangleq {\beta _{{\text{BR}}}}{{\text{e}}^{ - {\text{j}}2\pi {f_{\text{c}}}\tau _u^{{\text{BR}}}}}$为等效路径增益，${\beta _{{\text{BR}}}}$为RIS-BS信道的路径增益，$\tau _u^{{\text{BR}}} \triangleq \tau _{u,0,0}^{{\text{BR}}}$，$\vartheta$和$\omega$分别为RIS的离开角(Angle-Of-Departure, AOD)的仰角和方位角，$\psi$为BS的AOD，${{\boldsymbol{a}}_{{\text{BS}}}}\left( {\psi ,{f_q}} \right)$为BS的阵列导向矢量矩阵。${{\boldsymbol{a}}_{{\text{BS}}}}\left( {\psi ,{f_q}} \right)$表示为^[13]

定义

其中，${{\boldsymbol{H}}_{{\text{URB,}}u,q}}$为第$u$个UE在第$q$个子载波处的级联信道。

4.   两阶段级联信道估计方案

针对RIS辅助的宽带多用户系统的两阶段信道估计方案的主要思想如下：在第1阶段中，选择一个典型用户，其传输了${T_1}$个时隙的正交导频序列，通过一种基于梯度下降法的迭代优化算法来估计AOA、信道增益和时延参数，得到典型用户的级联信道。在第2阶段中，其他用户将发送${T_{{\text{up}}}} - {T_1}$个时隙的正交导频序列，该阶段利用RIS信道的空间相关性和用户信道之间的相关性减少信道估计的训练开销。

4.1   第1阶段：典型用户级联信道估计

4.1.1   问题描述

在这一阶段，选择最接近BS的UE作为典型用户，记为UE 1($u = 1$)。该阶段只有典型用户在${T_1}$个时隙中向BS发送正交导频信号，则BS在第$q$个子载波处的接收信号为

其中，${g_l} \triangleq {\overline \alpha _{_{{\text{RU,}}l}}}{\overline \beta _{{\text{BR}}}}$为第$l$条路径的级联信道增益，${\upsilon _{1,l}} \triangleq \tau _1^{{\text{BR}}} + \tau _{1,l}^{{\text{RU}}}$为UE 1的第$l$条路径的级联信道时延， ${\boldsymbol{W}} \triangleq {\text{blkdiag}}\left( {{{\boldsymbol{w}}_1},{{\boldsymbol{w}}_2}, \cdots ,{{\boldsymbol{w}}_{{T_1}}}} \right)$, ${{\boldsymbol{A}}_{1,q}}\left( {\vartheta ,\omega ,{f_q}} \right) \triangleq \left[ {{\left[ {{\boldsymbol{a}}_{{\text{RIS}}}^{\text{H}}(\vartheta ,\omega ,{f_q}){\text{diag}} \left( {{{\boldsymbol{\varPhi}} _1}} \right)} \right]}^{\text{T}}} ,{{\left[ {{\boldsymbol{a}}_{{\text{RIS}}}^{\text{H}}(\vartheta ,\omega ,{f_q}){\text{diag}} \left( {{{\boldsymbol{\varPhi}} _2}} \right)} \right]}^{\text{T}}}, \cdots \right.$, $\left.{{\left[ {{\boldsymbol{a}}_{{\text{RIS}}}^{\text{H}}(\vartheta ,\omega ,{f_q}){\text{diag}}\left( {{{\boldsymbol{\varPhi}} _{{T_1}}}} \right)} \right]}^{\text{T}}} \right]^{\text{T}}$, ${\overline {\boldsymbol{n}} _{1,q}} \triangleq \left[ {{ {{{\widetilde {\boldsymbol{n}}}^{\text{T}}_{1,q,1}}} }}\right.$, $\left.{{{{{\widetilde {\boldsymbol{n}}}^{\text{T}}_{1,q,2}}}}}, \cdots ,{{{{{\widetilde {\boldsymbol{n}}}^{\text{T}}_{1,q,{T_1}}}} }} \right]^{\text{T}}$。

BS处$Q$个子载波的接收信号为

其中，$\overline {\boldsymbol{W}} \triangleq {\text{blkdiag}}\left( {{\boldsymbol{W}},{\boldsymbol{W}}, \cdots ,{\boldsymbol{W}}} \right) \in {\mathbb{C}^{Q{T_1}N \times Q{T_1}N}}$, ${{\boldsymbol{N}}_1} \triangleq {\left[ {{{ {{{\overline {\boldsymbol{n}} }^{\text{T}}_{1,1}}} }} ,\; {{ {{{\overline {\boldsymbol{n}} }^{\text{T}}_{1,2}}} }} , \;\cdots , \;{{ {{{\overline {\boldsymbol{n}} }^{\text{T}}_{1,Q}}} }}} \right]^{\text{T}}}$, ${{\boldsymbol{\varLambda}} _l} \triangleq \Bigr[ \left[ {{\boldsymbol{A}}_{1,1}}(\vartheta ,\omega ,{f_1}) {{\boldsymbol{a}}_{{\text{RIS}}}}({\phi _l},{\theta _l},{f_1}) \otimes {{\boldsymbol{a}}_{{\text{BS}}}}(\psi ,{f_1}){{\text{e}}^{ - {\text{j}}2\pi {f_1}{\upsilon _{1,l}}}} \right]^{\text{T}}$, $\left[ {{\boldsymbol{A}}_{1,2}}(\vartheta ,\omega ,{f_2}) {{\boldsymbol{a}}_{{\text{RIS}}}}({\phi _l},{\theta _l},{f_2}) \otimes {{\boldsymbol{a}}_{{\text{BS}}}}(\psi , {f_2}){{\text{e}}^{ - {\text{j}}2\pi {f_2}{\upsilon _{1,l}}}} \right]^{\text{T}} , \cdots$, $\left.\left[ {{\boldsymbol{A}}_{1,Q}} (\vartheta ,\omega , {f_Q}) {{\boldsymbol{a}}_{{\text{RIS}}}}({\phi _l},{\theta _l},{f_Q}) \otimes {{\boldsymbol{a}}_{{\text{BS}}}}(\psi ,{f_Q}){{\text{e}}^{ - {\text{j}}2\pi {f_Q}{\upsilon _{1,l}}}} \right]^{\text{T}} \right]^{\text{T}}$。

由式(10)可知，级联信道的CSI由3种参数决定，即路径增益，角度信息和时延。此外，BS和RIS的位置固定，即$\vartheta$,$\omega$和$\psi$完全可以由其几何位置确定^[17]。因此，可将级联信道估计问题转化为参数恢复问题，即估计参数集$\left\{ {{g_l},{\phi _l},{\theta _l},{\upsilon _{1,l}}} \right\}_{l = 1}^L$。将式(10)改写为

其中，${\boldsymbol{g}} \triangleq {\left[ {{g_{_1}},{g_{_2}}, \cdots ,{g_{_{\text{L}}}}} \right]^{\text{T}}}$为信道增益矩阵，${\boldsymbol{K}} \triangleq \left[ {{\varLambda _1},{\varLambda _2}, \cdots ,{\varLambda _L}} \right]$为角度参数和时延参数组成的矩阵。

由于信道参数$\left\{ {{g_{_l}},{\phi _l},{\theta _l},{\upsilon _{1,l}}} \right\}_{l = 1}^L$的数量远小于${{\boldsymbol{Y}}_1}$的维数，可将式(11)的信道估计问题表示为

其中，$\widehat {\boldsymbol{g}} \triangleq {\left[ {{g_1},{g_1}, \cdots ,{g_{\widehat L}}} \right]^{\rm T}}$为估计的级联信道增益，$\widehat {\boldsymbol{\phi}} \triangleq {\left[ {{\phi _1},{\phi _2}, \cdots ,{\phi _{\widehat L}}} \right]^{\rm T}}$, $\widehat {\boldsymbol{\theta}} \triangleq {\left[ {{\theta _1},{\theta _2}, \cdots ,{\theta _{\widehat L}}} \right]^{\rm T}}$和$\widehat {\boldsymbol{\nu}} \triangleq \Bigr[ {v_{1,1}}, {v_{1,2}}, \cdots ,{v_{1,\widehat L}} \Bigr]^{\rm T}$为估计的角度参数矩阵，$\widehat L$为估计的路径数，$\varepsilon > 0$为误差容错参数。

4.1.2   基于梯度下降的迭代加权算法

求解式(12)的主要困难在于优化是一个NP-hard问题。文献[18]提出对数和函数对于稀疏信号恢复问题具有更强的稀疏性。利用对数和函数替换式(12)中的0范数，得到对数和最小化目标函数为

其中，$\xi > 0$是一个较小的正则化常数。

显然，式(13)中的优化问题是非凸的^[18]。通过添加数据拟合项将问题转化为如式(14)的无约束优化问题^[19]

其中，$\chi > 0$为控制稀疏性和数据拟合误差之间平衡的正则化参数。

为得到不增的目标函数值，利用优化最小化(Majorization Minimization, MM)方法求解式(14)，该方法通过迭代最小化一个简单的凸代理函数以优化给定的目标函数^[20]。借鉴文献[18]的思想，将$\varUpsilon \left( {{\boldsymbol{g}},{\boldsymbol{\phi}} ,{\boldsymbol{\theta}} ,{\boldsymbol{\nu}} } \right)$的代理函数建立为

其中，${\widehat {\boldsymbol{g}}^{\left( i \right)}}$为${\boldsymbol{g}}$第$i$次迭代时的估计值。

证明　为证明式(15)中$J\left( {{\boldsymbol{g}}|{{\widehat {\boldsymbol{g}}}^{\left( i \right)}}} \right)$为凸函数，令常数项${{C}} = \left| {{{\left[ g \right]}_l}} \right|_{\text{F}}^2 + \xi > 0$，变量$x = \left| {{{\left[ {{g^{\left( i \right)}}} \right]}_l}} \right|_{\text{F}}^2 + \xi > 0$，则有$k\left( x \right) = \displaystyle\sum\nolimits_{l = 1}^L {\left[ {\lg x + {{C}}/x - 1} \right]}$, ${k}'\left(x\right)= \left(x-{C}\right)/{x}^{3}$, ${k}''(x)=\left(2{C}-x\right)/{x}^{3}$，易证${k}''(x)\ge 0$，故$J\left( {{\boldsymbol{g}}|{{\widehat {\boldsymbol{g}}}^{\left( i \right)}}} \right)$为凸函数。证毕

将式(15)代入式(14)中，$\varUpsilon \left( {{\boldsymbol{g}},{\boldsymbol{\phi}} ,{\boldsymbol{\theta}} ,{\boldsymbol{\nu}} } \right)$的最小化等价于代理函数的最小化，即

其中，${G^{\left( i \right)}}\left( {{\boldsymbol{g}},{\boldsymbol{\phi}} ,{\boldsymbol{\theta}} ,{\boldsymbol{\nu}} } \right)$为第$i$次迭代时的迭代代理函数，${{\boldsymbol{D}}^{\left( i \right)}} = {\text{diag}}\left( 1\Big/\left( {{{\left| {\widehat g_{1,1}^{\left( i \right)}} \right|}^2} + \xi } \right),1\Big/\left( {{{\left| {\widehat g_{1,2}^{\left( i \right)}} \right|}^2} + \xi } \right) \cdots , 1\Big/\left( {{{\left| {\widehat g_{1,L}^{\left( i \right)}} \right|}^2} + \xi } \right) \right)$ 。

与文献[19]相似，可以优化式(16)中的${\boldsymbol{g}}$，以找到最佳点$\widehat {\boldsymbol{g}}$和对应的${G^{\left( i \right)}}$的最佳值，则有

将式(17)代入式(16)中，最优目标函数可表示为

最终，得到以最小化式(18)中的目标函数为目标的迭代算法。在所提参数估计算法中，在第$i$次迭代时，利用梯度下降法找到新的参数${\widehat {\boldsymbol{\phi}} ^{\left( {i + 1} \right)}}$, ${\widehat {\boldsymbol{\theta}} ^{\left( {i + 1} \right)}}$和${\widehat {\boldsymbol{\nu}} ^{\left( {i + 1} \right)}}$使得目标函数$G_{{\text{opt}}}^{\left( i \right)}\left( {{\boldsymbol{\phi}} ,{\boldsymbol{\theta}} ,{\boldsymbol{\nu}} } \right)$变得更小。在迭代过程中，动态正则化参数$\chi$随着迭代的进行不断更新，$\chi$更新方式为

其中，${\varTheta ^{\left( i \right)}}$为上一步的平方残差，$\zeta$为常数缩放因子，${\varTheta ^{\left( i \right)}} = \left\| {{{\boldsymbol{Y}}_1} - {{\overline {\boldsymbol{W}} }^{\text{H}}}{\boldsymbol{K}}\left( {{{\boldsymbol{\phi}} ^{\left( i \right)}},{{\boldsymbol{\theta}} ^{\left( i \right)}},{{\boldsymbol{\nu}} ^{\left( i \right)}}} \right){{\boldsymbol{g}}^{\left( i \right)}}} \right\|_{\text{F}}^2$。

考虑到实际场景中，$L$路径数是未知的。故所提信道估计算法中，将路径数$L$初始化为一个较大的值，且与真实路径数保持相同的数量级。在迭代过程中，将路径增益过小的路径视为噪声。基于此的一种剪枝策略为：如果$\left\| {g_l^{\left( {i + 1} \right)}} \right\|_{\text{F}}^2 < {g_{{\text{th}}}}, \forall l \in \left\{ {1,2, \cdots ,{{\widehat L}^{\left( i \right)}}} \right\}$, ${g_{{\text{th}}}}$为修剪阈值，则从${\widehat {\boldsymbol{g}}^{\left( {i + 1} \right)}}$中删除$g_{\widehat l}^{\left( {i + 1} \right)}$，从${\widehat {\boldsymbol{\phi}} ^{\left( {i + 1} \right)}}$, ${\widehat {\boldsymbol{\theta}} ^{\left( {i + 1} \right)}}$和${\widehat {\boldsymbol{\nu}} ^{\left( {i + 1} \right)}}$中删除$\widehat \phi _{\widehat l}^{\left( {i + 1} \right)}$, $\widehat \theta _{\widehat l}^{\left( {i + 1} \right)}$和$\widehat v_{1,\widehat l}^{\left( {i + 1} \right)}$，并更新${\widehat L^{\left( {i + 1} \right)}}$使其与${\widehat {\boldsymbol{\phi}} ^{\left( {i + 1} \right)}}$的长度相等。

基于梯度下降的迭代加权算法的流程如下：

步骤1　设置终止阈值${\varsigma _{{\text{th}}}}$，修剪阈值${g_{{\text{th}}}}$，根据式(17)初始化${\widehat {\boldsymbol{g}}^{\left( 0 \right)}} = {\boldsymbol{g}}_{{\text{opt}}}^{\left( 0 \right)}\left( {{{\widehat {\boldsymbol{\phi}} }^{\left( 0 \right)}},{{\widehat {\boldsymbol{\theta}} }^{\left( 0 \right)}},{{\widehat {\boldsymbol{\nu}} }^{\left( 0 \right)}}} \right)$。

步骤2　根据式(19)更新$\chi$。

步骤3　根据式(18)构造函数$G_{{\text{opt}}}^{\left( i \right)}\left( {{\boldsymbol{\phi}} ,{\boldsymbol{\theta}} ,{\boldsymbol{\nu}} } \right)$。

步骤4　利用梯度下降算法更新参数 ${\widehat {\boldsymbol{\phi}} ^{\left( {i + 1} \right)}}$, ${\widehat {\boldsymbol{\theta}} ^{\left( {i + 1} \right)}}$和 ${\widehat {\boldsymbol{\nu}} ^{\left( {i + 1} \right)}}$。

步骤5　根据式(17)估计路径增益${\widehat {\boldsymbol{g}}^{\left( {i + 1} \right)}}$。

步骤6　如果$\left\| {g_l^{\left( {i + 1} \right)}} \right\|_{\text{F}}^2 < {g_{{\text{th}}}}$，剪枝并更新路径数${\widehat L^{\left( {i + 1} \right)}}$。

步骤7　如果${\left\| {{{\widehat {\boldsymbol{g}}}^{\left( {i + 1} \right)}} - {{\widehat {\boldsymbol{g}}}^{\left( i \right)}}} \right\|_2} < {\varsigma _{{\text{th}}}}$，则返回步骤2；否则，得到$\widehat {\boldsymbol{g}}$, $\widehat {\boldsymbol{\phi}}$, $\widehat {\boldsymbol{\theta}}$, $\widehat {\boldsymbol{\nu}}$，由式(8)得到${\widehat {\boldsymbol{H}}_{{\text{URB,}}1,q}}, \forall q \in \left\{ {1,2, \cdots ,Q} \right\}$。

4.2   第2阶段：其他用户级联信道估计

由于每个RIS反射元件通过相同的RIS-BS信道将来自不同用户的信号反射到BS，所以可将任意用户的反射信道视为其他用户的反射信道的缩放版本^[21]。因此，在信道估计期间，只需要估计缩放因子(标量)，而不需要估计整个信道矢量。在该阶段，令${H_{{\text{RU}},u,q,m}}$表示第$u$个UE到RIS的第$m$个反射元件在第$q$个子载波处的频域信道响应，${\mu _{u,m}}$表示与RIS的第$m$个反射元件相关的缩放因子，利用不同UE的UE-RIS信道之间的强相关性^[14]，可得

对整个${{\boldsymbol{H}}_{{\text{RU}},u,q}}$信道，则有

其中，${{\boldsymbol{\mu}} _{u,q}} \triangleq {\left[ {{\mu _{u,1,q}},{\mu _{u,2,q}}, \cdots ,{\mu _{u,M,q}}} \right]^{\text{T}}}$为与${{\boldsymbol{H}}_{{\text{RU}},u,q}}$相关的缩放矢量。

由式(8)和式(21)可得

其中，${{\boldsymbol{H}}_{{\text{URB,}}1,b,q}} \triangleq {{\boldsymbol{H}}_{{\text{URB,}}1,q}}{\text{diag}}\left( {{{\boldsymbol{\varPhi}} _b}} \right)$。

利用式(21)和式(22)的关系，$\left\{ {{{\boldsymbol{H}}_{{\text{URB,}}u,q}}} \right\}_{u = 2}^U$可以由${{\boldsymbol{H}}_{{\text{URB,}}1,q}}$和$\left\{ {{{\boldsymbol{\mu}} _{u,q}}} \right\}_{u = 2}^U$的信息估计得到。因此，通过估计$\left\{ {{{\boldsymbol{H}}_{{\text{URB,}}1,q}},\left\{ {{{\boldsymbol{\mu}} _{u,q}}} \right\}_{u = 2}^U} \right\}$而不是估计$\left\{ {{{\boldsymbol{H}}_{{\text{URB,}}u,q}}} \right\}_{u = 1}^U$，从而显著减少待估计的未知信道系数的数量。将式(22)代入式(1)并向量化，在$U - 1$个用户发送了${T_2} = {T_{{\text{up}}}} - {T_1}$时隙的正交导频序列后，BS接收到第$q$个子载波的信号为

其中，${{\boldsymbol{\varXi}} _q} \triangleq \left[ {{\left[ {{{\boldsymbol{B}}_{2,1,q}},{{\boldsymbol{B}}_{2,2,q}}, \cdots ,{{\boldsymbol{B}}_{2,{T_2},q}}} \right]}^{\text{T}}}, \cdots , \left[ {{\boldsymbol{B}}_{U,1,q}}, {{\boldsymbol{B}}_{U,2,q}}, \cdots , {{\boldsymbol{B}}_{U,{T_2},q}} \right]^{\text{T}} \right]$, ${{\boldsymbol{B}}_{u,b,q}} \triangleq {\boldsymbol{s}}_{u,b,q}^{\text{T}} \otimes \left( {{\boldsymbol{w}}_b^{\text{H}}{{\boldsymbol{H}}_{{\text{URB,}}1,b,q}}} \right)$, ${{\boldsymbol{\varGamma}} _q} \triangleq {\left[ {{{\boldsymbol{\mu}} _{2,q}},{{\boldsymbol{\mu}} _{3,q}}, \cdots ,{{\boldsymbol{\mu}} _{U,q}}} \right]^{\text{T}}}$, ${{\boldsymbol{z}}_q} \triangleq \left[ {\text{vec}}\left( {{{\boldsymbol{n}}_{1,q}}} \right),{\text{vec}}\left( {{{\boldsymbol{n}}_{2,q}}} \right), \cdots , {\text{vec}}\left( {{{\boldsymbol{n}}_{{T_2},q}}} \right) \right]^{\text{T}}$。

BS处$Q$个子载波的接收信号为

其中，${\boldsymbol{\varpi}} \triangleq {\text{blkdiag}} \left[ {{{\boldsymbol{\varXi}} _1},{{\boldsymbol{\varXi}} _2}, \cdots ,{{\boldsymbol{\varXi}} _Q}} \right] \in {\mathbb{C}^{{T_2}PQ{N_{{\text{RF}}}} \times \left( {U - 1} \right)QM}}$, ${\boldsymbol{\varPsi}} \triangleq {\left[ {{{\boldsymbol{\varGamma}} _1},{{\boldsymbol{\varGamma}} _2}, \cdots ,{{\boldsymbol{\varGamma}} _P}} \right]^{\text{T}}} \in {\mathbb{C}^{\left( {U - 1} \right)QM \times 1}}$, ${\boldsymbol{Z}} \triangleq \left[ {{\boldsymbol{z}}_1},{{\boldsymbol{z}}_2}, \cdots , {{\boldsymbol{z}}_Q} \right]^{\text{T}}$。

文献[22]表明，使得式(24)存在唯一解的条件为：(1)${T_2}PQ{N_{{\text{RF}}}} \ge \left( {U - 1} \right)QM$；(2) ${\text{rank}}\left( {\boldsymbol{\varpi}} \right) = \left( {U - 1} \right)QM$。因此，该阶段传输的导频数应满足：${T_2}P \ge \left\lceil {\left( {U - 1} \right)M/{N_{{\text{RF}}}}} \right\rceil$。由LS算法可得

其中，$\widehat {\boldsymbol{\varPsi}}$为${\boldsymbol{\varPsi}}$的估计值。

5.   最小导频开销与计算复杂度分析

5.1   最小导频开销分析

在第1阶段中，基于梯度下降的迭代加权算法估计了$L$个级联信道增益$\widehat {\boldsymbol{g}}$，及$3L$个的角度参数$\widehat \phi$,$\widehat \theta$和$\widehat \nu$，因此，在该阶段的导频开销为${T_1} \ge 4L$；在第2阶段中，因为需要确保式(24)存在唯一解，所需最小导频数目为$\left\lceil {\left( {U - 1} \right)M/{N_{{\text{RF}}}}} \right\rceil$。因此，提出的两阶段级联信道估计方案最小导频开销为$4LP + \left\lceil {\left( {U - 1} \right)M/{N_{{\text{RF}}}}} \right\rceil$。

对比方案中，文献[4]最小导频开销为$UL\lg M$，文献[5]最小导频开销为$UL\lg \left( {{L^2}{E_{{\text{UE}}}}{E_{{\text{RIS}}}}} \right)$，其中${E_{{\text{UE}}}}$和${E_{{\text{RIS}}}}$为角分辨率，文献[6]最小导频开销为$U{T_{{\text{up}}}}P\left\lceil {N/{N_{{\text{RF}}}}} \right\rceil \left\lceil {M/{M_{{\text{RF}}}}} \right\rceil$，其中${M_{{\text{RF}}}}$为RIS处RF个数，文献[8]最小导频开销为$M + \max \left( U - 1, \left\lceil {\left( {U - 1} \right)M/{N_{{\text{RF}}}}} \right\rceil \right)$，文献[14]最小导频开销为$2\left( {M + 1} \right)/\rho + U\left\lceil {M/N} \right\rceil + U$，其中，$\rho \gg 1$表示一个常量系数。相比之下，所提方案的最小导频开销小于文献[4]和文献[6]。

5.2   计算复杂度分析

在第1阶段的基于梯度下降的迭代加权算法中，复杂度主要来自步骤4中梯度计算，计算复杂度为$O\left( {Q{T_1}{N_{{\text{RF}}}}{L^2}\left( {Q{T_1}N + L} \right)} \right)$；第2阶段利用了LS算法，复杂度主要来自式(25)中的伪逆计算，计算复杂度为$O\left( {Q^3}{U^3}{M^3} + 2{Q^3}{U^2}{M^2}{T_2}P{N_{{\text{RF}}}} + {Q^2}{T_2}P{N_{{\text{RF}}}}UM \right)$。因此，所提两阶段级联信道估计方案总复杂度为$O\left( {Q{T_1}{N_{{\text{RF}}}}{L^2}\left( {Q{T_1}N + L} \right)} \right) + O\left( {Q^3}{U^3}{M^3} + 2{Q^3}{U^2}{M^2}{T_2}P{N_{{\text{RF}}}} + {Q^2}{T_2}P{N_{{\text{RF}}}}UM\right)$。

对比方案中，文献[4]主要计算复杂度为$O\left( {{t^3}U \;+\; 2{P^3}{t^2}U{T_2}P{N_{{\text{RF}}}} \;+\; {Q^2}U{T_2}P{N_{{\text{RF}}}}t} \right)$+$O\left( {{Q^2}U{T_2}P{N_{{\text{RF}}}}t} \right) + O\left( {{Q^2}U{T_2}P{N_{{\text{RF}}}}I} \right)$，其中$t$为迭代次数，$I$为量化网格数，文献[5]主要计算复杂度为$O\left( {2{t^3}UQTP + 4{Q^2}{t^2}{T^2}{P^2}UN + 2{Q^2}{T^2}{P^2}UNt} \right)$ + $O\left( {{T^2}{P^2}QUN{E_{{\text{UE}}}}{E_{{\text{RIS}}}}{E_{{\text{BS}}}}} \right)$ + $O\left( {{T^2}{P^2}QU{E_{{\text{BS}}}}} \right)$，其中${E_{{\text{BS}}}}$为角分辨率，文献[6]主要计算复杂度为$O\left( {UQTP{{\left( {M + N} \right)}^{3.5}}} \right)$，文献[8]主要计算复杂度为$O \left( {{Q^3}{N^3}{M^3} + 2{Q^3}{N^2}{M^2}{T_1}P{N_{{\text{RF}}}} + {P^2}NM{T_1}P{N_{{\text{RF}}}}} \right)$+ $O\left( {{Q^3}{U^3}{M^3} + 2{Q^3}{U^2}{M^2}{T_2}P{N_{{\text{RF}}}} + {Q^2}{T_2}P{N_{{\text{RF}}}}UM} \right)$，文献[14]主要计算复杂度为$O\left( {Q^3}{N^3}{M^3}U + 2{Q^3}{N^2}{M^2}U{T_2}P{N_{{\text{RF}}}} + {Q^2}U{T_2}P{N_{{\text{RF}}}}NM \right)$ + $O\left( {M^3} + {M^2}NL + MNLt \right)$。虽然本文提出的两阶段级联信道估计方案在第2阶段也涉及了伪逆计算，但在第1阶段中，因为路径数$L$较小，所以复杂度$O\left( {P{T_1}{N_{{\text{RF}}}}{L^2}\left( {Q{T_1}N + L} \right)} \right)$较小。因此，相比于对比方案，所提方案的总的复杂度更低。

6.   仿真分析

本节探究所提两阶段级联信道估计方案的可行性和优越性。参数设置如下：BS处的天线数$N = 256$，RIS处的天线数$M = 256$，UE数量$U = 8$，UE和RIS之间的路径数$L = 3$，传输带宽$B = 20{\text{ GHz}}$，载波频率${f_{\text{c}}} = 150{\text{ GHz}}$，子载波数$Q = 256$，信噪比(Signal to Noise Ratio, SNR)为–5 dB。与本文所提两阶段信道估计方案的对比算法中，选择未考虑波束斜视影响的RIS辅助通信的窄带和宽带信道估计方案：基于压缩感知的分布式正交匹配追踪算法的文献[4]、基于压缩感知的3阶段信道估计方案的文献[5]、基于混合RIS结构的2阶段信道估计方案的文献[6]、多用户信道估计方案的文献[8]和基于坐标下降信道估计的文献[14]。

图1给出了不同SNR下不同信道估计方案的NMSE性能对比。由于文献[4]和文献[5]中的信道估计方案从离散网格中选择估计角度，网格失配使得其性能不如其他算法，且网格字典大小的限制使得随着SNR的增加估计精度逐渐达到饱和。相比之下，所提2阶段级联信道估计方案中基于梯度下降的迭代加权算法可以避免量化误差，NMSE性能明显好于对比的5种信道估计方案。

图  1  不同信噪比下的归一化均方误差

下载: 全尺寸图片 幻灯片

图2给出了不同信道估计方案的NMSE性能与传输带宽$B$的关系。由图2可知，随着传输带宽$B$的增加，波束斜视的影响逐渐增大，文献[4]和文献[5]中的信道估计方案的性能越来越差。这是因为，对于较大的传输带宽$B$，不同子载波处的宽带波束空间信道的支撑集差异性更大，导致公共稀疏信道支撑集假设不再成立。文献[6]、文献[8]和文献[14]中的信道估计方案的性能不随传输带宽$B$的增加而降低，因为这3种信道估计方案是逐个估计每个子载波对应的信道。相反，所提2阶段级联信道估计方案对传输带宽$B$是鲁棒的。

图  2  不同传输带宽下的归一化均方误差

下载: 全尺寸图片 幻灯片

图3给出了不同UE数量$U$下各种信道估计方案的信道估计性能关系。由图3可知，随着UE数量的逐渐增加，所提信道估计方案的性能有所下降，当UE数量大于20时，所提方案的性能不如文献[6]和文献[14]的信道估计方案。这是因为，所提2阶段级联信道估计方案的性能与典型用户的估计准确度密切相关，随着用户数量的增加，其他用户利用信道间的强相关性进行信道估计的信道误差将会逐渐增加。

图  3  不同UE数量下的归一化均方误差

下载: 全尺寸图片 幻灯片

图4给出了不同信道估计方案的导频长度$P{T_{{\text{up}}}}$与NMSE性能的关系。由图4可知，随着导频符号长度的增加，所有信道估计方案的准确度都会提高，因为通过大量的测量值可以获得更好的恢复和估计精度。其次，对于相同导频长度$P{T_{{\text{up}}}}$，所提两阶段级联信道估计方案始终比其他信道估计方案具有更好的NMSE性能。

图  4  不同导频数量下的归一化均方误差

下载: 全尺寸图片 幻灯片

7.   结束语

本文研究了波束斜视影响下RIS辅助宽带多用户通信系统中的上行信道估计问题，并利用THz信道的稀疏性，提出一种两阶段级联信道估计方案。首先，将信道估计问题转化为无约束优化问题，并利用基于梯度下降法的迭代优化算法来估计AOA、信道增益和时延参数；然后，基于信道间的强相关性估计其他用户的级联反射信道。仿真结果表明，与传统宽带算法相比，本算法具有更好、更稳定的信道估计性能。针对信道估计中用户数增多导致性能下降和矩阵伪逆的计算导致较高复杂度的问题，后续工作将展开对用户进行分组和降低复杂度的分析与研究。