OTFS信道估计原理与MATLAB实现:从时延-多普勒域到嵌入式导频
2026/9/13 14:49:27 网站建设 项目流程

简介:针对正交时频空间(OTFS)在大规模MIMO系统下的信道估计难题,这份MATLAB仿真资源为无线通信领域的研究者与工程师提供了一套可直接运行的完整工程。程序涵盖OTFS与OFDM符号生成、大规模MIMO信道建模、迫零与最小均方误差检测、基于正交匹配追踪的稀疏信道估计等核心算法,并支持单发单收、多发单收等仿真场景,能够灵活评估不同信噪比、移动速度下的误码率与归一化均方误差,绘制多组性能对比曲线,便于深入验证算法效果。压缩包共68个文件,以61个Matlab脚本为主,覆盖参数配置、信道生成、收发链路、检测与绘图全流程;另有3个txt说明文档、2个C语言辅助源文件、1个fig结果图与1个mat数据文件,包体大小为31.39MB。txt文件用于算法说明与使用指引,C文件用于C-MEX接口加速,fig和mat保存了典型性能曲线与仿真数据,整体结构清晰,便于按模块学习或二次开发。目前已有1185人学习下载,适合正在从事OTFS信道估计科研实验、毕业设计或课程设计的学生作为实战参考。

1. 高速移动下OFDM变脆,OTFS的信道估计也得换个打法

OFDM在LTE和Wi-Fi里几乎统治了无线物理层,但把场景换成高铁、低轨卫星这类高速移动环境,子载波间隔和信道时间变化之间的冲突就会把系统性能拉下来。OTFS把调制符号直接放在时延-多普勒域上,让每个符号经历完整的时频选择性信道,理论上能把这部分损失压回去。可它给信道估计出的难题是:时延和多普勒都是连续变量,需要同时估计出二维冲激响应,而不是像OFDM那样在每个子载波上放一串导频就完事。配合大规模MIMO之后,天线数一上来,导频开销和估计复杂度又被放大一个数量级。下面面向已经在MATLAB里跑过链路仿真、想真正把OTFS信道估计模块落地的人,从DD域模型、导频设计到可跑的估计代码,把参数设法和常见坑一次说清。

2. 时延-多普勒域输入输出关系与嵌入式导频的约束条件

2.1 OTFS发送链路的域变换:从DD域符号到时域波形

OTFS的发送端可以拆成两步理解。第一步,把原始调制符号排列成一个M行N列的二维网格X_dd,行方向叫时延采样,列方向叫多普勒采样。这个网格就是DD域。第二步,对这个网格做逆辛傅里叶变换(ISFFT),映射到传统的时频域X_tf,再按OFDM的方式做Heisenberg变换生成时域波形。

ISFFT的具体方向在MATLAB里最容易搞反。我常用的约定是:沿时延维做M点FFT,沿多普勒维做N点IFFT。

X_tf = fft(X_dd, [], 1) / sqrt(M); % 时延维FFT X_tf = ifft(X_tf, [], 2) * sqrt(N); % 多普勒维IFFT

这里fft(...,[],1)的作用对象是矩阵的每一列,对应时延维;ifft(...,[],2)作用在每一行,对应多普勒维。缩放因子/sqrt(M)*sqrt(N)是为了让能量在两次变换间守恒,实际实现时只要收发两端对称,缩放怎么分配不影响估计结果。

接收端做逆过程:先Wigner变换得到Y_tf,再做SFFT回到DD域。方向与发送端完全相反,即先沿时延维IFFT,再多普勒维FFT。

2.2 双选择性信道的时延-多普勒稀疏表示

无线信道的时延-多普勒响应可以写成P条可分辨路径的叠加:

h(τ,ν) = Σ h_i δ(τ-τ_i) δ(ν-ν_i)

每条路径有一个复增益h_i、时延τ_i、多普勒ν_i。低速场景下ν_i接近0,OFDM只需要对抗频率选择性;高速场景下ν_i变大,OFDM子载波间的正交性被破坏,性能急剧下滑。OTFS的基函数本身就是时延-多普勒域上的脉冲,所以信道在这个域里是稀疏且慢变的,每个抽头只落在少数几个格点上。

离散化时,时延索引和多普勒索引由分辨率和物理量共同决定:

  • 时延分辨率 = 1/(M·Δf),对应索引 l_i = round(τ_i · MΔf)
  • 多普勒分辨率 = 1/(N·T),对应索引 k_i = round(ν_i · NT)

其中Δf是子载波间隔,T是一个OTFS符号的持续时间(不含CP)。这两个分辨率直接决定信道估计能不能分辨出相邻路径。

2.3 嵌入式导频、保护带尺寸与导频开销

OTFS信道估计最常用的方案是嵌入式导频:在DD网格上放一个功率经过放大的导频符号,周围留出一块空白的保护带。保护带的作用是防止数据符号的二维循环卷积能量泄漏进导频检测窗口。

保护带尺寸由最大时延索引L_max和最大多普勒索引K_max决定。时延维度上需要2L_max+1个格点,多普勒维度上需要2K_max+1个格点。导频开销的粗略公式是:

开销 = (2L_max+1)(2K_max+1) / (M·N)

参数符号OTFS典型值受谁约束
时延维度M128~1024带宽与CP长度
多普勒维度N16~128多普勒分辨率
最大时延索引L_max4~32多径时延扩展
最大多普勒索引K_max2~8多普勒扩展
导频功率放大β²10~20 dB峰值检测阈值与噪声底

导频功率放大β²是OTFS估计里一个不起眼但决定成败的参数。放大倍数太小,峰值会被噪声淹没;放大倍数太大,会挤占发射机的功率预算,降低数据符号的有效信噪比。我一般取13~17 dB,再根据仿真里峰值旁瓣的情况微调。

2.4 为什么OFDM的频域导频方案不能直接搬过来

OFDM的导频在时频格点上只被当前格点的信道响应H[n,m]相乘,估计的是每个资源格上的复增益。OTFS不是这样,接收端Y_dd是发送X_dd与信道冲激响应的二维循环卷积,一个导频符号的能量会扩散到所有时延-多普勒抽头上。如果照搬OFDM在每个格点放导频的思路,导频之间会互相干扰,保护带也失去意义。所以OTFS必须把导频放在DD域,并在二维平面上做峰值搜索。

3. MATLAB实现OTFS帧结构:从ISFFT方向到二维峰值估计

3.1 帧参数与DD域发送网格构造

先把仿真参数定义清楚,后面所有代码都依赖这一组变量。

% ---- OTFS帧参数 ---- N = 32; % 多普勒域符号数 M = 128; % 时延域采样数 cpLen = 16; % 循环前缀长度(样本数),必须大于最大时延 fs = M * 15e3; % 采样率:15kHz作为子载波间隔 T = M / fs; % 一个OTFS符号持续时间(不含CP) % ---- 导频与保护带参数 ---- Lwin = 2; % 时延维保护半径 Kwin = 2; % 多普勒维保护半径 lp = 20; kp = 6; % 导频位置(时延索引,多普勒索引) beta = sqrt(20); % 导频幅度放大,约13 dB % ---- 生成DD域发送网格 ---- rng(42); dataSym = qammod(randi([0 3], M*(N-1), 1), 4, 'gray'); X_dd = zeros(M, N); X_dd(:, 2:N) = reshape(dataSym, M, N-1); % 数据放在多普勒维第2~N列 % ---- 嵌入导频并清空保护带 ---- X_dd(lp-Lwin:lp+Lwin, kp-Kwin:kp+Kwin) = 0; X_dd(lp, kp) = beta;

数据占满第2到第N列,导频周围的5×5窗口被显式清零,再放入放大后的导频符号。这样保护带内没有数据符号,二维循环卷积后不会污染导频位置。这里kp=6Kwin=2,保护带不会越界;如果kp太小,保护带会循环绕到DD网格另一端,仿真结果会莫名出现多余峰值。

3.2 发送端:ISFFT与Heisenberg变换

% ---- ISFFT:时延维FFT,多普勒维IFFT ---- X_tf = fft(X_dd, [], 1) / sqrt(M); X_tf = ifft(X_tf, [], 2) * sqrt(N); % ---- Heisenberg变换:逐列OFDM调制并串接 ---- tx = zeros((M + cpLen) * N, 1); for n = 1:N sym = ifft(X_tf(:, n), M); % OFDM符号调制 blk = [sym(end-cpLen+1:end); sym]; % 加循环前缀 tx((n-1)*(M+cpLen)+1 : n*(M+cpLen)) = blk; end tx = tx / sqrt(mean(abs(tx).^2)); % 发送功率归一化

逐列做IFFT生成时域符号,再拼接成整帧。这里用循环而不是矩阵化,是为了让CP插入逻辑更直观,M和N规模不大时性能足够。功率归一化放在最后统一做,避免导频放大影响整体功率统计。注意这里ISFFT顺序必须与接收端SFFT严格对称,否则解调出来的是噪声。

3.3 时变信道与时延-多普勒抽头实现

% ---- 信道抽头参数 ---- P = 4; tauIdx = [0 3 8 14]; % 时延,单位:采样周期 nuIdx = [0 -1 1 2]; % 多普勒,单位:多普勒分辨率 alpha = [0.85 0.60 0.45 0.30]; % 路径幅度 phi = [0.3 -0.7 1.1 -0.4]; % 路径相位 vmax = 1 / (N * T); % 多普勒分辨率 f_d = nuIdx * vmax; % 实际多普勒频偏 % ---- 时域逐径叠加 ---- ht = zeros(size(tx)); for i = 1:P delay = tauIdx(i); dphase = exp(1j*2*pi*f_d(i)*(0:length(tx)-1)'/fs); shifted = [zeros(delay, 1); tx(1:end-delay)]; ht = ht + alpha(i)*exp(1j*phi(i)) * shifted .* dphase; end % ---- 加噪声 ---- snr_dB = 15; noisePwr = 10^(-snr_dB/10); rx = ht + sqrt(noisePwr/2) * (randn(size(ht)) + 1j*randn(size(ht)));

每径的时延通过移位实现,多普勒通过逐样本相位旋转实现。相位旋转作用在整个OTFS帧上,等效于对每条路径叠加一个连续频偏。这块的采样率fs必须和发射端一致,否则多普勒频偏的离散相位会算错。

3.4 接收端解调与二维峰值信道估计

% ---- Wigner变换:逐OFDM块解调 ---- Y_tf = zeros(M, N); for n = 1:N blk = rx((n-1)*(M+cpLen)+1 : n*(M+cpLen)); Y_tf(:, n) = fft(blk, M); end % ---- SFFT:先时延维IFFT,再多普勒维FFT ---- Y_dd = ifft(Y_tf, [], 1) * sqrt(M); Y_dd = fft(Y_dd, [], 2) / sqrt(N); % ---- 在导频窗口内做峰值检测 ---- winL = lp-Lwin:lp+Lwin; winK = kp-Kwin:kp+Kwin; PeakBox = Y_dd(winL, winK); th = 0.3 * max(abs(PeakBox(:))); % 相对峰值阈值 mask = abs(PeakBox) > th; [rows, cols] = find(mask); l_est = rows - 1 + winL(1); % 还原DD域绝对索引 k_est = cols - 1 + winK(1); h_est = PeakBox(mask) / beta; % 除以导频幅度得信道增益

find返回的rows是时延维坐标,cols是多普勒维坐标,还原后得到的就是每条路径在DD网格上的时延索引和多普勒索引。h_est是复信道增益,去除了导频放大的影响。整个估计过程没有用任何遍历或搜索,核心只有一次阈值比较。

3.5 参数速查:从物理量到MATLAB变量

物理量公式MATLAB变量本例取值
时延分辨率1/(MΔf)1/fs0.52 μs
多普勒分辨率1/(NT)vmax468.75 Hz
最大时延索引τ_max·MΔfmax(tauIdx)14
最大多普勒索引ν_max·NTmax(abs(nuIdx))2
保护带尺寸(2L_max+1)(2K_max+1)(2Lwin+1)(2Kwin+1)25

4. 大规模MIMO天线阵下的导频复用与低复杂度估计结构

4.1 大规模MIMO下信道估计的维度结构

大规模MIMO把单天线OTFS链路扩展成多天线版本后,待估计对象从M×N的二维矩阵变成A×M×N的三维张量,A是基站天线数。上行场景中,每个用户发送同一帧OTFS信号,基站每根天线独立做SFFT,得到各自的Y_dd(:,:,a)。训练开销由用户的导频决定,不会因为天线数增加而线性增长,这一点和OFDM大规模MIMO的结论一致。

4.2 多用户嵌入式导频复用规则

多用户场景下,导频复用最常见的方式是给不同用户分配不同的时延-多普勒位置。用户u的导频位置要满足两个条件:自身保护带完整,并且与其他用户的保护带不重叠。

U = 4; tauMax = max(tauIdx); nuMax = max(abs(nuIdx)); lp_u = zeros(1, U); kp_u = zeros(1, U); for u = 1:U lp_u(u) = tauMax + 2 + (u-1) * (2*tauMax + 4); kp_u(u) = nuMax + 2; end

时延维上每个用户错开(2·tauMax+4)个格点,多普勒维共用同一位置。这样每个用户导频周围都有完整的矩形保护带,路径时延在tauMax以内时不会跨用户泄漏。实际系统中还应该考虑路径多普勒超过nuMax的扩散,因此多普勒维的间隔常见做法是留出2·(nuMax+1)。

4.3 逐天线估计的复杂度控制

估计复杂度主要由峰值搜索的扫描范围决定。按3.4节的阈值检测,每根天线只需要在U个导频窗口内做比较,复杂度是O(U·(2Lwin+1)(2Kwin+1)),与M·N无关。这比全网格最小二乘或消息传递低得多,也是嵌入式导频在大规模MIMO里受欢迎的根本原因。

当天线数A=64、用户数U=4、窗口25个格点时,总检测点数只有6400,MATLAB里一次矩阵索引和比较就完成。如果改用基于压缩感知的OMP或AMP重构信道,复杂度会随路径数P和迭代次数上升,但好处是可以在导频保护带更小的情况下工作,适合导频开销敏感的场景。天线越多,时延-多普勒支撑越一致,这类稀疏重构算法的增益越明显。

4.4 从估计结果到预编码用的完整信道

峰值检测得到的是每条路径的时延、多普勒和复增益,但预编码需要完整的时频域信道H[n,m,a]。常见做法是把估计结果重新投射回时频域:

H_est_TF = zeros(M, N, A); for u = 1:U for path = 1:nPeak lt = l_est(path); kt = k_est(path); H_est_TF(:, :, a) = H_est_TF(:, :, a) + ... h_est(path) * exp(1j*2*pi*kt*(0:N-1)/N) * ... % 多普勒维 exp(-1j*2*pi*lt*(0:M-1).'/M).'; % 时延维 end end

这段构造用的是时延-多普勒二维正弦基叠加,物理含义是每条路径对每个时频格点都贡献一个相位旋转。实际预编码前还应该做插值和去噪,但对链路级评估来说,这个投射已经可以支撑波束成形增益的计算。

5. 分数多普勒、阈值选择与OTFS信道估计的验证方法

5.1 用归一化MSE验证估计精度

验证估计器好坏,我习惯先算归一化MSE:

H_true = zeros(M, N); for i = 1:P lt = mod(lp + tauIdx(i) - 1, M) + 1; kt = mod(kp + nuIdx(i) - 1, N) + 1; H_true(lt, kt) = alpha(i) * exp(1j*phi(i)); end % 稀疏向量化后比较 nMSE = sum(abs(h_est(:) - H_true(mask)).^2) / sum(abs(H_true(mask)).^2);

注意h_est的长度是检测到的路径数,如果阈值漏检,MSE会显著变大。单次仿真的MSE波动很大,至少跑50次取平均才可信。

5.2 分数多普勒泄漏:过采样SFFT

真实信道的多普勒频偏很少正好落在整数倍分辨率上。当nuIdx=0.5时,能量会泄漏到相邻多普勒bin,峰值旁瓣变高,阈值检测容易把一条径误判成多条。一个有效的处理办法是在SFFT前对多普勒维做零填充过采样:

Kos = 4; % 过采样倍数 Y_tf_pad = [Y_tf, zeros(M, N*(Kos-1))]; Y_dd_pad = ifft(Y_tf_pad, [], 1) * sqrt(M); Y_dd_pad = fft(Y_dd_pad, [], 2) / sqrt(Kos*N);

过采样把多普勒网格细化到1/(Kos·N·T),泄漏能量重新集中到同一个bin上。代价是Y_dd维度变成M×(Kos·N),峰值扫描的计算量也跟着乘Kos。我一般只在验证阶段开过采样,批量仿真时保持整数多普勒。

5.3 自适应阈值的快速检查法

固定相对阈值th=0.3·max()在信噪比变化时不够稳。一个简单实用的自检方案是利用保护带空白区域估计噪声底:

noiseBox = Y_dd(lp-Lwin-6:lp+Lwin+6, kp-Kwin-6:kp+Kwin+6); noiseBox(ismember(...)) = []; % 剔除导频窗口部分 noiseFloor = mean(abs(noiseBox(:))); th = noiseFloor + 3 * std(abs(noiseBox(:)));

这个阈值相当于噪声底之上3个标准差,信噪比在5到25 dB范围内都能保持较低的漏检率。如果路径幅度差异超过20 dB,弱径会被阈值滤掉,此时可以把阈值下调到2个标准差,再通过最小间距合并相邻峰值。验证时如果看到估计出的路径数大于P,先查是不是保护带尺寸不够,再查阈值和过采样,三个地方都调过之后基本能收敛到正确结果。

本文还有配套的精品资源,点击获取

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询