1.   引 言

可重构智能表面(Reconfigurable Intelligent Surface, RIS)作为第6代(Sixth Generation, 6G)移动通信潜在技术之一，具有广阔的应用前景^[1]。RIS辅助大规模多输入多输出(Multiple Input Multiple Output, MIMO)系统可以显著提高无线通信覆盖率、吞吐量和频谱效率^[2]。毫米波信道的高路损极易阻断用户(User Equipment, UE)与基站(Base Station, BS)的视距(Line of Sight, LoS)径，严重影响通信质量。引入RIS动态调整信号幅度或相移可以提高传输性能^[3]。文献[4]利用无人机辅助RIS提高其反射能力及覆盖范围，可以最大化中继网络传输速率。文献[5]在信源与用户间放置RIS可以降低保密中断概率，增强系统保密性能。实现上述目标需要获取UE-RIS-BS链路的准确信道状态信息(Channel State Information, CSI)。

在RIS辅助的毫米波通信系统中，一般BS和RIS的分布是固定的，其链路为准静态信道，多个时隙内仅需估计1次。而环境突然变化(如车辆和行人堵塞)以及终端缓慢变化会导致UE-RIS信道突变和角度变化。为保证估计性能，需要在多个时隙内追踪UE-RIS信道角度^[6]并检测信道突变。目前，RIS分为无源和有源两种架构，文献[7]针对无源RIS架构提出最小二乘(Least Square, LS)算法以估计上行链路级联信道，但是该架构在RIS端无法获得观测数据，所以信道估计和追踪只能在BS或UE端完成，既降低了RIS波束成形增益，又具有相当大的计算复杂度。文献[8]为降低计算复杂度，采用高斯赛德尔(Gauss Seidel, GS)迭代算法进行估计，但在UE角度变化时估计灵活性较差。文献[9]提出一种基于稀疏矩阵分解的信道估计方法，但需要对每个RIS无源元件进行单独的幅度控制，导致成本极为昂贵。文献[10]在时变场景下将级联信道建模为状态空间模型以捕获信道相关性，并采用卡尔曼滤波(Kalman Filter, KF)算法估计时变信道。文献[11]使用KF分别追踪LoS信道和级联信道，进一步优化导频和RIS反射系数矩阵以提高估计精度。RIS配置射频(Radio Frequency, RF)链组成有源RIS架构，可以在RIS端处理信号并分离估计UE-RIS-BS信道，利用无线信道固有稀疏性能以少量训练开销提取特定信道参数，降低计算复杂度^[12,13]。如文献[14]针对有源RIS结构，利用波束空间域信道稀疏性，提出一种多目标交替优化信道估计方案，可以在RIS侧灵活估计信道并节省训练开销，但不涉及波束空间域角度变化，不利于实际系统的应用。文献[15]将RIS连接RF链用于基带测量，利用深度学习方法估计RIS信道，证明了信道感知时间对最优RIS相移矩阵的影响，但是需要大量数据集进行长时间离线训练。文献[16]利用有限RF链分别估计UE-RIS信道和RIS-BS信道，提出一种基于求根多重信号分类(Root-MUltiple SIgnal Classification, Root-MUSIC)算法的信道估计方案，可以获得级联链路独立CSI，有效提升估计灵活性，但是Root-MUSIC算法涉及求解多项式的根，仍然存在较大计算复杂度。

针对RIS辅助的毫米波通信系统，本文提出一种基于牛顿算法的信道追踪方案。所提方案将RIS连接“L”型RF链，首先使用2维快速傅里叶变换(Two-Dimensional Fast Fourier Transform, 2D-FFT)算法分离估计第1个时隙UE-RIS和RIS-BS信道角度参数^[17]，并使用最大似然(Maximum Likelihood, ML)算法分离估计级联信道路径增益完成初始化，最后在后续时隙使用牛顿算法追踪UE-RIS链路角度，并对UE-RIS信道进行突变检测，如果检测到突变则再次初始化，否则继续使用牛顿算法追踪角度。仿真结果验证该方案在环境突然变化和UE缓慢移动时具有优越性能，并且分离估计级联信道可以减少估计开销，经过分析，所提方案计算复杂度远小于LS算法、线性最小均方误差(Linear Minimum Mean Squared Error, LMMSE)算法和KF算法^[18]，极大节约算力资源，实现计算复杂度与性能之间的良好折中。

2.   系统模型

本文研究RIS辅助上行无线通信系统毫米波信道追踪问题，如图1所示，BS配备${N_{\text{R}}}$根均匀线性阵列天线(Uniform Linear Array, ULA)，服务配备${N_{\text{T}}}$根天线的单用户。RIS具有${N_{\text{I}}} = N_{\text{I}}^{\text{h}} \times N_{\text{I}}^{\text{v}}$个反射单元，$N_{\text{I}}^{\text{h}}$和$N_{\text{I}}^{\text{v}}$为RIS水平和竖直方向反射单元数目。由于RIS无法进行信号处理，本文以RIS板右下角反射单元为原点，向水平与竖直方向分别延伸连接$L$个RF链构成“L型”分布。为降低系统能耗，RF链数远小于RIS板反射单元数($L < < N_{\text{I}}^{\text{h}}$,$L < < N_{\text{I}}^{\text{v}}$)，且UE-BS的LoS径被阻断，仅能通过RIS辅助完成通信。

图 1  RIS辅助上行无线通信系统

下载: 全尺寸图片 幻灯片

该系统BS和RIS空间位置保持不变，RIS-BS信道${\boldsymbol{G}}$为准静态信道，当终端缓慢变化时会导致UE-RIS信道${\boldsymbol{H}}(t)$突变和角度变化，所以将UE-RIS信道${\boldsymbol{H}}(t) \in {\mathbb{C}^{{N_{\text{I}}} \times {N_{\text{T}}}}}$和RIS-BS信道${\boldsymbol{G}} \in {\mathbb{C}^{{N_{\text{R}}} \times {N_{\text{I}}}}}$分别建模为

其中，${L_a}$和${L_b}$为${\boldsymbol{H}}(t)$和${\boldsymbol{G}}$路径数，${a_{{l_a}}}(t)$为${\boldsymbol{H}}(t)$第${l_a}$条路径增益，${b_{{l_b}}}$为${\boldsymbol{G}}$第${l_b}$条路径增益，${(·)}^{\text{H}}$表示共轭转置，$\varphi _{\boldsymbol{H}}^{{l_a}}(t)$为UE离开角(Angle of Departure, AoD)，$\theta _{\boldsymbol{H}}^{{l_a}}(t)$和$\phi _{\boldsymbol{H}}^{{l_a}}(t)$为RIS到达角(Angle of Arrival, AoA)中仰角和方位角，${\gamma _{{l_b}}}$和${\varphi _{{l_b}}}$为RIS的AoD中仰角和方位角，${\psi _{{l_b}}}$为BS的AoA，${\boldsymbol{H}}(t)$导向矢量${\boldsymbol{a}}_{\text{U}}^{}(\varphi _{\boldsymbol{H}}^{{l_a}}(t))$和${\boldsymbol{a}}_{{\text{UR}}}^{}(\theta _{\boldsymbol{H}}^{{l_a}}(t),\phi _{\boldsymbol{H}}^{{l_a}}(t))$为

其中，${(\cdot)}^{\text{T}}$表示转置，定义${u_{{l_a}}}(t)$和${v_{{l_a}}}(t)$及${w_{{l_a}}}(t)$分别为

并定义${d_{\text{U}}}$和${d_{\text{I}}}$为 UE和RIS的天线间距，$\lambda$为波长，导向矢量${{\boldsymbol{a}}_{\text{B}}}({\psi _{{l_b}}})$与式(3)类似，将$\varphi _{\boldsymbol{H}}^{{l_a}}(t)$和${d_{\text{U}}}$用${\psi _{{l_b}}}$和${d_{\text{B}}}$代替，${d_{\text{B}}}$为BS天线间距，${{\boldsymbol{a}}_{{\text{RB}}}}({\gamma _{{l_b}}},{\varphi _{{l_b}}})$与式(4)类似，将$\theta _{\boldsymbol{H}}^{{l_a}}(t)$和$\phi _{\boldsymbol{H}}^{{l_a}}(t)$用${\gamma _{{l_b}}}$和${\varphi _{{l_b}}}$代替。RIS未连接RF链的水平导向矢量为${\boldsymbol{a}}_{{\text{UR}}}^{\text{h}}(\phi _{\boldsymbol{H}}^{{l_a}}(t)) = [1, \cdots ,{{\text{e}}^{{\text{j2}}\pi m{v_{{l_a}}}(t)}}, \cdots ,{{\text{e}}^{{\text{j}}2\pi (N_{\text{I}}^{\text{h}} - 1){v_{{l_a}}}(t)}}]^{\text{T}}$，RIS未连接RF链的竖直导向矢量为${\boldsymbol{a}}_{{\text{UR}}}^{\text{v}}(\theta _{\boldsymbol{H}}^{{l_a}}(t)) = [1, \cdots ,{{\text{e}}^{{\text{j2}}\pi n{w_{{l_a}}}(t))}}, \cdots , {{\text{e}}^{{\text{j}}2\pi (N_{\text{I}}^{\text{v}} - 1){w_{{l_a}}}(t))}}]^{\text{T}}$，BS在第$t$个时隙的接收信号${\boldsymbol{Y}}(t) \in {\mathbb{C}^{{N_{\text{R}}} \times T}}$为

其中，${\boldsymbol{\varGamma}} (t) \in {\mathbb{C}^{{N_{\text{T}}} \times T}}$为UE发送的导频信号，满足正交条件${\boldsymbol{\varGamma}} (t){\boldsymbol{\varGamma}}^{\text{H}} (t) = {{\boldsymbol{I}}_{{N_{\text{T}}}}}$，$T$为时隙数，${\boldsymbol{N}}(t) \in {\mathbb{C}^{{N_{\text{R}}} \times T}} \sim \mathcal{C}\mathcal{N}{\text{(}}0,\sigma _k^2{\text{)}}$是均值为0，方差为$\sigma _k^2$的复高斯白噪声。在不同时隙信号传输过程，${\boldsymbol{G}}$为准静态信道，仅需估计1次，而${\boldsymbol{H}}(t)$可能发生角度变化和信道突变，需要频繁追踪角度并检测突变以保证性能。

第$t$个时隙UE-RIS所有路径角度参数集合${\boldsymbol{O}}(t)$为

相邻时隙的角度变化可建模为1阶马尔可夫链^[18]，所以第$t$个时隙存在角度变化的角度模型为

其中，${\boldsymbol{u}}(t){\rm{～}}\mathcal{C}\mathcal{N}(0,\sigma _u^2)$是均值为0，方差为$\sigma _u^2$的复高斯白噪声，并定义${{\boldsymbol{Q}}_{\text{U}}}$为${\boldsymbol{u}}(t)$的协方差矩阵。

3.   信道追踪

由式(1)和式(2)说明获得角度和路径增益可恢复整个信道矩阵。本文在第1个时隙与${\boldsymbol{H}}(t)$突变的时隙使用2D-FFT初始化估计角度，后续时隙利用牛顿算法追踪${\boldsymbol{H}}(t)$角度，并使用ML估计每个时隙路径增益。

3.1   初始化信道参数

本节利用2D-FFT在第1个时隙与${\boldsymbol{H}}(t)$发生突变的时隙得到信号空间谱，估计UE-RIS和RIS-BS信道角度，然后利用ML算法估计路径增益。

3.1.1   初始化角度参数

因为RIS在水平和竖直方向上估计角度的方法一致，本文选择以RIS所连RF链水平方向接收信号为例来分析角度参数初始化。RIS水平方向接收信号${\boldsymbol{Y}}_{\text{R}}^{\text{h}}(t) \in {\mathbb{C}^{L \times T}}$为

其中，${\boldsymbol{H}}^{\text{h}}(t) \in {\mathbb{C}^{L \times {N_{\text{T}}}}}$为UE到RIS水平方向信道矩阵，${\boldsymbol{N}}^{\text{h}}(t) \in {\mathbb{C}^{L \times T}}$为信号经过${\boldsymbol{H}}^{\text{h}}(t)$叠加的高斯白噪声，将接收信号与导频解相关后得到${\boldsymbol{Y}}_{\text{R}}^{\text{h}'}(t) \in {{\mathbb{C}}^{L \times {N_{\text{T}}}}}$为

其中，${\boldsymbol{N}}_{}^{{\text{h}}'} (t) = {\boldsymbol{N}}^{\text{h}}(t){\boldsymbol{\varGamma}}^{\text{H}} {(t)}$为等效噪声，对${\boldsymbol{Y}}_{\text{R}}^{{\text{h}}'}(t)$补0得到${\boldsymbol{Y}}_{\text{R}}^{{\text{h}}' } (t) \in {\mathbb{C}^{L' \times {N'_{\text{T}}}}}$，$L'$和${N'_{\text{T}}}$为补0后行与列的数目，再通过2D-FFT得到信号空间谱$\tilde {\boldsymbol{Y}}_{\text{R}}^{\text{h}}(t) \in {\mathbb{C}^{L' \times {N'_{\text{T}}}}}$为

其中，$\tilde {\boldsymbol{H}} ^{\text{h}}{{(t)}} \in {\mathbb{C}^{L' \times {N'_{\text{T}}}}}$为${\boldsymbol{H}}_{}^{\text{h}}(t)$的空间谱，${\tilde {\boldsymbol{N}}^{\text{h}}}(t) \in {\mathbb{C}^{L' \times {N'_{\text{T}}}}}$为${\boldsymbol{N}}^{{\text{h}}'}(t)$的空间谱，${{\boldsymbol{F}}_{{{L'}_{}}}} \in {\mathbb{C}^{L' \times L'}}$和${\boldsymbol{F}}_{{N'_{\text{T}}}}^{} \in {\mathbb{C}^{{N'_{\text{T}}} \times {N'_{\text{T}}}}}$为补0后傅里叶矩阵，$\tilde {\boldsymbol{H}}_{}^{\text{h}}{{(t)}}$第$m$行$n$列的元素$\tilde {\boldsymbol{H}}_{m,n}^{\text{h}}(t)$为

其中，定义$f(x) = ({1 \mathord{\left/ {\vphantom {1 {M)}}} \right. } {M)}} \cdot ({{\sin (\pi xM)} \mathord{\left/ {\vphantom {{\sin (\pi xM)} {\sin (\pi x)}}} \right. } {\sin (\pi x)}}) \cdot {{\exp( - {\rm{j}}}}2\pi x({{M - 1} \mathord{\left/ {\vphantom {{M - 1} 2}} \right. } 2}))$，经过谱峰搜索确定第${l_a}$条路径峰值点坐标为$(m_{{l_a}}^{\text{h}}(t),n_{{l_a}}^{\text{h}}(t))$，根据式(14)可得估计值${\hat u_{{l_a}}}(t) = {{m_{{l_a}}^{\text{h}}(t)} \mathord{\left/ {\vphantom {{m_{{l_a}}^{\text{h}}(t)} {L'}}} \right. } {L'}}$, ${\hat v_{{l_a}}}(t) = {{n_{{l_a}}^{\text{h}}(t)} \mathord{\left/ {\vphantom {{n_{{l_a}}^{\text{h}}(t)} {{N'_{\text{T}}}}}} \right. } {{N'_{\text{T}}}}}$，所以UE-RIS水平方向第${l_a}$条路径角度估计值为

将UE-RIS所连RF链竖直方向接收信号进行2D-FFT可估计出$\hat \theta _{\boldsymbol{H}}^{{l_a}}(t)$为

其中，$n_{{l_a}}^{\text{v}}(t)$为$\hat \theta _{\boldsymbol{H}}^{{l_a}}(t)$对应峰值点坐标，同理初始化信道${\boldsymbol{G}}$的角度参数。

3.1.2   初始化路径增益参数

使用ML算法初始化路径增益，将式(11)改写为

其中，${\boldsymbol{H}}^{\text{h}}(t) = {\boldsymbol{A}}_{{\text{UR}}}^{\text{h}}(t){\boldsymbol{B}}(t){{\boldsymbol{A}}_{\text{U}}^{\text{H}}(t)}$，定义${\boldsymbol{H}}_{}^{\text{h}}(t)$阵列导向矩阵${\boldsymbol{A}}_{{\text{UR}}}^{\text{h}}(t) \in {\mathbb{C}^{L \times {L_a}}}$, ${\boldsymbol{A}}_{\text{U}}^{}(t) \in {\mathbb{C}^{{N_{\text{T}}} \times {L_a}}}$和路径增益矩阵${\boldsymbol{B}}(t) \in {\mathbb{C}^{{L_a} \times {L_a}}}$分别为

利用ML算法求解第$t$个时隙UE-RIS信道矩阵路径增益估计值${{\hat {\boldsymbol{B}}}}{(t)_{{\text{ML}}}}$为^[16]

其中，${(\cdot)}^{-1}$为求逆运算，$\hat {\boldsymbol{A}}_{{\text{UR}}}^{\text{h}}(t)$和$\hat {\boldsymbol{A}}_{\text{U}}^{}(t)$为${\boldsymbol{A}}_{{\text{UR}}}^{\text{h}}(t)$和${\boldsymbol{A}}_{\text{U}}^{}(t)$的估计值，同理初始化信道${\boldsymbol{G}}$路径增益参数。

3.2   追踪$H(t)$角度参数

本节利用牛顿算法^[19,20]迭代追踪后续时隙${\boldsymbol{H}}(t)$的角度。牛顿算法每次迭代追踪1条路径的所有角度，然后将该路径分量从带噪的叠加信号中移除得到残差，最后从残差中继续追踪下个路径的角度参数。第$t$个时隙第${l_a} - 1$次迭代结束后的残差${\boldsymbol{R}}$为

其中，$i = 0,1,\cdots,{l_a} - 1$, ${l_a} = 1,2,\cdots,{L_a}$, $\{ {\hat a_i}(t),{\hat u_i}(t), {\hat v_i}(t)\} _{i = 0,1,\cdots,{l_a} - 1}$为前${l_a} - 1$次迭代所得到的估计值，再进行第${l_a}$次迭代得到新路径角度参数估计值，需使残差能量$E({u_{{l_a}}}(t),{v_{{l_a}}}(t))$最小化，即

其中，$|| \cdot ||_{\text{F}}^{}$为矩阵${\text{F}}$范数，最小化$E({u_{{l_a}}}(t),{v_{{l_a}}}(t))$需将其泰勒展开为

其中，定义${\hat u_{{l_a}}}(t - 1)$为${u_{{l_a}}}(t)$在第$t - 1$个时隙的估计值，${\hat v_{{l_a}}}(t - 1)$为${v_{{l_a}}}(t)$在第$t - 1$个时隙的估计值。并定义${{\dot {\boldsymbol{E}}}}({\hat u_{{l_a}}}(t - 1),{\hat v_{{l_a}}}(t - 1))$和${{\ddot {\boldsymbol{E}}}}({\hat u_{{l_a}}}(t - 1), {\hat v_{{l_a}}}(t - 1))$分别为

其中，${E'_{{{\hat u}_{{l_a}}}(t - 1)}}({\hat u_{{l_a}}}(t - 1),\;{\hat v_{{l_a}}}(t - 1))$为$E({\hat u_{{l_a}}}(t - 1), \;{\hat v_{{l_a}}}(t - 1))$对${\hat u_{{l_a}}}(t - 1)$求1阶偏导，${E''_{{{\hat u}_{{l_a}}}(t - 1) {{\hat u}_{{l_a}}}(t - 1)}} ({\hat u_{{l_a}}}(t - 1),{\hat v_{{l_a}}}(t - 1))$为$E({\hat u_{{l_a}}}(t - 1),{\hat v_{{l_a}}}(t - 1))$对${\hat u_{{l_a}}}(t - 1)$求2阶偏导，${E''_{{{\hat u}_{{l_a}}}(t - 1){{\hat v}_{{l_a}}}(t - 1)}}({\hat u_{{l_a}}}(t - 1), {\hat v_{{l_a}}}(t - 1))$为$E({\hat u_{{l_a}}}(t - 1),{\hat v_{{l_a}}}(t - 1))$先对${\hat u_{{l_a}}}(t - 1)$求偏导后再对${\hat v_{{l_a}}}(t - 1)$求偏导，其余变量定义同理。对式(25)左右两边同时求偏导可得${E'_{{u_{{l_a}}}(t)}}({u_{{l_a}}}(t), {v_{{l_a}}}(t))$和${E'_{{v_{{l_a}}}(t)}}({u_{{l_a}}}(t),{v_{{l_a}}}(t))$分别为

分别令式(28)和式(29)为0即可得到$({u_{{l_a}}}(t),{v_{{l_a}}}(t))$的取值，但是第${l_a}$次追踪过程仅追踪1次得到的角度值不精确，需要多次循环。假设循环${R_{\text{c}}}$次，则第${l_a}$次追踪第$k$次循环估计的${\hat u_{{l_a}}}{(t)^{(k)}}$和${\hat v_{{l_a}}}{(t)^{(k)}}$分别为

其中，$k = 1,2, \cdots ,{R_{\text{c}}}$, ${\hat u_{{l_a}}}{(t)^{(0)}} = {\hat u_{{l_a}}}(t - 1)$, ${\hat v_{{l_a}}}{(t)^{(0)}} = {\hat v_{{l_a}}}(t - 1)$。利用式(30)循环迭代追踪角度，然后使用式(22)估计路径增益，恢复出第$t$个时隙的${\boldsymbol{H}}^{\text{h}}(t)$。同理恢复出UE-RIS竖直方向的信道矩阵。

3.3   检测信道${\boldsymbol{H(t) }}$突变

本节对信道矩阵${\boldsymbol{H}}(t)$是否突变进行检测。存在两种假设，${{{H}}_{\text{1}}}$：信道矩阵突变，${{{H}}_{\text{0}}}$：信道矩阵未突变。定义似然函数为

其中，${{\hat {\boldsymbol{H}}}}_{}^{\text{h}}(t) = \displaystyle\sum\nolimits_{{l_a} = 1}^{{L_a}}{\text{ }} {\hat a_{{l_a}}}(t){\boldsymbol{a}}_{{\text{UR}}}^{\text{h}}({\hat v_{{l_a}}}(t)){\boldsymbol{a}}_{\text{U}}^{\text{H}}( {{\hat u_{{l_a}}}(t)}$为信道${\boldsymbol{H}}^{\text{h}}(t)$的估计值，${\boldsymbol{Q}}_{\text{U}}^{}$为噪声协方差矩阵，${\det ^{{N_{\text{T}}}}}({{\boldsymbol{Q}}_{\text{U}}})$为计算${{\boldsymbol{Q}}_{\text{U}}}$行列式的${N_{\text{T}}}$次方，当发生${H_{\text{0}}}$事件时，式(31)${\boldsymbol{Y}}_{\text{R}}^{{\text{h}}'}(t) - \hat {\boldsymbol{H}}_{}^{\text{h}}(t)$将趋近于${\boldsymbol{0}}$，根据纽曼-皮尔逊准则^[21]得到对数似然比为

当判决为${H_{\text{1}}}$事件需满足

其中，$\gamma$是预先定义的判决门限，预设发生${H_{\text{0}}}$事件被判决为${H_{\text{1}}}$事件的虚警概率${P_{{\text{fa}}}}$为

给定虚警概率${P_{{\text{fa}}}}$可确定门限$\gamma$的取值，因为$2\ln (L({\boldsymbol{Y}}_{\text{R}}^{{\text{h}}'} (t)))$服从卡方分布$\chi _{2L{N_{\text{T}}}}^2$，定义${Q_{\chi _v^2}}(x) = \displaystyle\int_x^\infty {p(t){\text{d}}t}$为右尾概率，其中$p(t)$为概率密度函数，可得门限$\gamma$为

依据以上判决准则，当判决为${H_{\text{1}}}$事件可认为信道突变，此时重新利用2D-FFT初始化信道矩阵，开始下一次追踪。综上，本文信道追踪方案算法流程如算法1所示。

算法1 基于牛顿算法的信道追踪方案
输入：时隙数：$t = 1,2,\cdots,I$；导频信号：${\boldsymbol{\varGamma}} (t)$；RIS接收信号：${\boldsymbol{Y}}_{\text{R} }^{\text{h} }(t)$,${\boldsymbol{Y}}_{\text{R} }^{\text{v} }(t)$；BS接收信号：${\boldsymbol{Y}}(t)$；循环次数：${R_{\text{c}}}$
输出：UE-RIS信道估计值${{\hat {\boldsymbol H}}}(t)$；RIS-BS信道估计值${{\hat {\boldsymbol G}}}$
for ${l_a} = 1,2,\cdots,{L_a}$
(1)${\text{when }}t = 1$，2D-FFT求解角度：根据式(15)、式(16)、式(17)分别计算$\hat \varphi _{\boldsymbol{H}}^{{l_a}}(1),\hat \phi _{\boldsymbol{H}}^{{l_a}}(1),\hat \theta _{\boldsymbol{H}}^{{l_a}}(1)$，同理求解${\psi _{{l_b}}},{\gamma _{{l_b}}},{\varphi _{{l_b}}}$ 　　　　　　　　 ML求解路径增益：根据式(22)计算${a_{{l_a}}}(1)$，同理求解${b_{{l_b}}}$
(2)for $t = 2,3,\cdots,I$
for $k = 1,2,\cdots,{R_{\text{c}}}$
$\left[ \begin{gathered} {{\hat u}_{{l_a}}}{(t)^{(k)}} \\ {{\hat v}_{{l_a}}}{(t)^{(k)}} \\ \end{gathered} \right] = \left[ \begin{gathered} {{\hat u}_{{l_a}}}{(t)^{(k - 1)}} \\ {{\hat v}_{{l_a}}}{(t)^{(k - 1)}} \\ \end{gathered} \right] - {{\ddot {\boldsymbol E}}}({\hat u_{{l_a}}}(t - 1),{\hat v_{{l_a}}}(t - 1)) \cdot {{\dot {\boldsymbol E}}}({\hat u_{{l_a}}}(t - 1),{\hat v_{{l_a}}}(t - 1))$
end for
${\hat u_{{l_a}}}(t) = {\hat u_{{l_a}}}{(t)^{({R_{\text{c}}})}},{\hat v_{{l_a}}}(t) = {\hat v_{{l_a}}}{(t)^{({R_{\text{c}}})}}$
求解角度：根据式(15)、式(16)分别计算$\hat \varphi _{\boldsymbol{H}}^{{l_a}}(t),\hat \phi _{\boldsymbol{H}}^{{l_a}}(t)$，同理求解$\hat \theta _{\boldsymbol{H}}^{{l_a}}(t)$
估计增益：根据式(22)计算${a_{{l_a}}}(t)$
信道突变检测：若$\ln (L({\boldsymbol{Y}}_{\text{R} }^{ {\text{h} }' }(t))) > \gamma$；返回步骤1
end for
end for 　(3)输出信道估计值$\hat {\boldsymbol{H} }(t) = \displaystyle\sum\nolimits_{ {l_a} = 1}^{ {L_a} } { { {\hat a}_{ {l_a} } }(t){\boldsymbol{a} }_{ {\text{UR} } }^{}(\hat \theta _{\boldsymbol{H} }^{ {l_a} }(t),\hat \phi _{\boldsymbol{H} }^{ {l_a} }(t)){ {{\boldsymbol{a} }_{\text{U} }^{\text{H} }(\hat \varphi _{\boldsymbol{H} }^{ {l_a} }(t))} } } {\text{, } }\hat {\boldsymbol{G} }{\text{ = } }\displaystyle\sum\nolimits_{ {l_b} = 1}^{ {L_b} } { { {\hat b}_{ {l_b} } }{ {\boldsymbol{a} }_{\text{B} } }({ {\hat \psi }_{ {l_b} } }){ {{\boldsymbol{a} }_{ {\text{RB} } }^{\text{H} }({ {\hat \gamma }_{ {l_b} } },{ {\hat \varphi }_{ {l_b} } })} } }$

下载: 导出CSV 

| 显示表格

4.   计算复杂度分析

本文信道追踪计算复杂度由追踪${\boldsymbol{H}}(t)$角度和估计路径增益两部分组成。作为对比，引入文献[18]的KF算法追踪本文场景信道${\boldsymbol{H}}(t)$角度，并在第1个时隙和${\boldsymbol{H}}(t)$突变的时隙使用2D-FFT初始化角度，利用ML估计路径增益，同时引入LS和LMMSE算法进行对比。本文将算法的乘法和加法运算量之和作为计算复杂度的衡量指标，将4种算法的复杂度对比结果列于表1中。

表 1  不同算法计算复杂度对比

所用算法	计算复杂度
LS算法	$O({N_{\text{T}}}^3 + N_{\text{T}}^2(2(T - 1) + 1 + 2N_{\text{I}}^{}) + 2N_{\text{I}}^{}N_{\text{T}}^{}(T - 1) + {N_{\text{T}}})$
LMMSE算法	$O(2N_{\text{T}}^3 + N_{\text{T}}^2(2(T - 1) + 1 + 8N_{\text{I}}^{}) + 2N_{\text{I}}^{}N_{\text{T}}^{}(T - 2) + 4N_{\text{T}}^{})$
KF算法	$O(16L_a^3 + {(L{N_{\text{T} } })^3} + {(L{N_{\text{T} } })^2}(4{L_a} + 1) + 2{L_a}^2(4L{N_{\text{T} } } + 3{N_{\text{T} } } + 2L + 1) + 2L{N_{\text{T} } }({L_a} + 1) - {L_a}L)$
本文算法	$O(8{L_a}^3 + 2{L_a}^2(3{N_{\text{T} } } + 2L - 1) + ({L_a} - 1)(2L{N_{\text{T} } } + {N_{\text{T} } }(4L - 2) + 22{R_{\text{c} } }) + L{N_{\text{T} } }(2T - 1) - {L_a}L) \\$

下载: 导出CSV 

| 显示表格

LS算法和LMMSE算法计算复杂度仅与RIS反射元件数$N_{\text{I}}^{}$、时隙数$T$及发送天线数${N_{\text{T}}}$相关，KF算法与本文算法还与RF链数$L$、路径数${L_a}$及循环次数${R_{\text{c}}}$相关。值得一提的是，KF算法涉及矩阵求逆以及高维矩阵相乘，因此其计算复杂度包含多个立方项$O(16L_a^3 + {(L{N_{\text{T}}})^3})$。将上述算法应用于RIS反射单元数$N_{\text{I}}^{} = 40$、水平/竖直方向上RF链数$L = 5$、循环次数${R_{\text{c}}} = 3$、用户天线数${N_{\text{T}}} = 5$、路径数${L_a} = 3$和时隙数$T = 10$的场景，KF算法复杂度最高为26460，是本文算法的17倍，LMMSE算法复杂度为11945，是本文算法的8倍，LS算法复杂度为6205，是本文算法的4倍。综上，验证了本文信道追踪方案可以极大减少计算复杂度。

5.   仿真结果与分析

本节通过仿真验证所提方案在RIS通信系统中的性能表现，仿真参数设置如表2所示。

表 2  主要仿真参数设置

参数	数值
UE端天线数${N_{\text{T}}}$	$5$
RIS端反射单元数${N_{\text{I}}}$	$40$
BS端天线数${N_{\text{R}}}$	$64$
RIS连接水平/竖直RF链数$L$	$5$
UE-RIS路径数${L_a}$	$1$
载波频率${f_{\text{c}}}$	$30{\text{ GHz}}$
UE端离开角$\hat \varphi _{\mathbf{H}}^{{l_a}}$	${50^ \circ }$
RIS端仰角$\theta _{\mathbf{H}}^{{l_a}}$和水平角$\phi _{\mathbf{H}}^{{l_a}}$	${30}^{\circ }和{20}^{\circ }$

下载: 导出CSV 

| 显示表格

采用归一化均方误差衡量所提信道追踪方案的性能，表达式为

其中，$|| \cdot ||_{\text{F}}^{}$为矩阵${\text{F}}$范数，${\boldsymbol{H}}$以${\boldsymbol{H}}(t)$或者${\boldsymbol{G}}$代替为真实信道，$\hat {\boldsymbol{H}}$以$\hat {\boldsymbol{H}}(t)$或者$\hat {\boldsymbol{G}}$代替为估计信道。

为更好地验证所提方案性能，引入文献[10]所用LS和LMMSE算法作为对比算法，因为该文献场景与本文类似，同时也引入文献[18]用于信道追踪的KF算法进行对比。图2为不同算法NMSE性能随SNR变化的对比曲线，设置角度变化方差$\sigma _u^2 = 0.5$。仿真结果表明，任意SNR下，本文方案性能均优于LS与LMMSE算法，因为本文方案将角度追踪转化为极值求解，若相邻时隙角度变化量较小，则使用上一时隙追踪的角度进一步优化当前时隙角度，虽相比于KF算法性能较弱，但计算复杂度远远小于KF算法，所以本文方案可以在角度变化场景以极低计算复杂度达到较好的追踪性能。

图 2  不同算法NMSE与SNR的对比曲线

下载: 全尺寸图片 幻灯片

图3为不同算法NMSE性能随角度变化方差$\sigma _u^2$变化的对比曲线，设置SNR=20 dB。仿真结果表明，本文方案与KF算法NMSE性能随$\sigma _u^2$增加而降低，因为它们均利用前一时隙信道参数相关性，$\sigma _u^2$越大对追踪性能影响越严重。而LS和LMMSE算法NMSE性能几乎保持不变，因为它们仅与当前时隙发送导频信号和接收信号相关，与角度变化方差无关。

图 3  不同算法NMSE与角度变化方差的对比曲线

下载: 全尺寸图片 幻灯片

图4为不同算法NMSE性能随时隙变化的对比曲线，设置SNR=20 dB, $\sigma _u^2 = 0.5$。仿真结果表明，当追踪时隙增加时，本文方案和KF算法NMSE性能逐渐降低，因为它们在本时隙追踪角度与上一时隙角度参数相关，会在多个时隙内产生误差传播，但在前3个时隙仍优于LS和LMMSE算法，验证了本文方案的可行性。而LS和LMMSE算法整体NMSE性能基本保持不变，因为它们仅利用当前时隙接收信号和导频。值得一提的是，本文方案计算复杂度在所有算法中最低且性能表现优越，可以很好地实现低复杂度信道追踪，有效节约计算资源，达到性能与复杂度的良好折中。

图 4  不同算法NMSE与时隙的对比曲线

下载: 全尺寸图片 幻灯片

6.   结束语

针对RIS辅助的毫米波上行无线通信系统，本文提出一种基于牛顿算法的低复杂度信道追踪方案。所提方案首先基于2D-FFT算法初始化信道角度，然后基于ML算法估计路径增益。在后续时隙基于牛顿算法追踪UE-RIS信道角度参数，并且对UE-RIS信道矩阵进行检测，判决其是否突变，若信道突变则再次初始化。本文通过将RIS部分反射元素连接RF链，分离估计UE-RIS和RIS-BS级联信道，可以有效降低估计开销。相比于其他算法，本文方案计算复杂度最低，能够极大节省算力资源，并且也未因计算复杂度降低而损失太多性能。仿真结果表明，本文方案具有优越的性能表现，可以在复杂度与性能之间取得很好的平衡。