JADE算法：角度与时延联合估计的MATLAB实现与避坑指南

简介：这是一套角度和时延联合估计（JADE）算法的MATLAB实现程序，面向无线通信、雷达探测、声学成像等方向的研究者与工程师，用于从多通道接收信号中同时估计到达角度和时延。程序以5个m脚本文件构成，压缩包约5KB，覆盖信号建模、二维FFT变换、四元数变换、独立成分分析（ICA）以及角度/时延提取等核心环节，代码结构紧凑，便于对照学习与二次修改。已有648人浏览学习。通过阅读与运行这套程序，读者可理解JADE从多通道信号模型到参数估计的完整流程，掌握ICA分离、复共轭对称性利用及优化求解等关键知识，也能将其应用于实际阵列信号处理任务，快速验证算法效果。

1. JADE 是什么：一次测量拿到每条路径的角度和到达时间

做阵列信号处理或信道测量时，经常遇到这样的问题：MUSIC 谱能给出几个来波方向，相关峰能给出几条时延，但角度和时延是分别估计出来的。室内多径、毫米波传播场景里，一条路径的到达角和时延往往同时影响接收数据，分开估计不仅精度受损，还会出现“角度对了但时延对应错路径”的配对难题。角度和时延联合估计（JADE）算法要解决的就是这个：从一次频域快照里同时估计出每条多径的到达角和传播时延，并且利用同一组子空间结构完成自动配对。这篇文章以一套可复现的 MATLAB 程序为主线，从频域模型、仿真信号生成，到二维 ESPRIT 实现、参数设定和踩坑记录，一路写到实测前需要留意的几个习惯。适合做信道探测、室内定位、测向和超分辨估计的工程师，也适合正在复现经典 JADE 算法但被配对步骤卡住的同学。

2. 先立住模型：两个旋转不变性如何把角度和时延“硬绑”在一起

JADE 的核心不在“估计”而在“联合”。如果每一步都独立估计，再做配对，误差会逐级传递；如果直接在二维参数空间搜索，计算量又太大。我常用的切入点是频域快照模型，它把时延变成频率维上的相位坡，把角度变成空间维上的相位坡，最后通过两个旋转不变性方程同时提取参数。

2.1 频率快照模型：时延如何表现为频域相位坡

考虑一个 M 元均匀线阵（ULA），阵元间距为 d，接收 L 条多径信号。发射端发送已知导频，接收端在每个阵元上做频域采样，得到 Nf 个频点的信道频域响应。第 l 条路径的到达角为 θ_l，时延为 τ_l，复增益为 β_l。频点 f_k 上的接收数据可以写为：

X(f_k) = Σ_{l=1}^{L} β_l a(θ_l) e^{-j2π f_k τ_l} + N(f_k)

其中 a(θ_l) 是 ULA 导向矢量。这个式子的关键在 e^{-j2π f_k τ_l}：时延不再是一个独立待估计的时间偏移，而是频率维度上的一条线性相位坡。坡的斜率就是 τ_l，所以只要测出频率维的相位变化率，就能反推出时延。

把 Nf 个频点的观测按列堆叠成一个 M×Nf 矩阵 X，然后向量化，得到每一条路径对应的联合导向矢量：

s_l = vec( β_l a(θ_l) e(τ_l)^T ) = β_l e(τ_l) ⊗ a(θ_l)

其中 e(τ_l) = [1, e^{-j2πΔf τ_l}, …, e^{-j2π(Nf-1)Δf τ_l}]^T。如果收集 K 个导频符号的频域快照，每个符号的复增益 β_l 不同，就构成一个 MNf×K 的数据矩阵 Y。这个模型的好处是：角度和时延被编码进同一个张量积结构里，后面提取参数时可以共享同一个信号子空间，配对天然就是联合的。

2.2 两个平移不变性：空间维旋转算子和频率维旋转算子

ULA 的导向矢量有一个经典性质：相邻阵元之间只差一个固定相位

μ_l = 2π (d/λ) sin θ_l

如果把联合导向矢量 s_l 在空间维上截掉最后一个阵元，和截掉第一个阵元相比，前者相差 e^{jμ_l}。这就是空间维的旋转不变性。

频域维同样存在这样的结构。等间隔频率采样下，联合导向矢量中相邻频点之间相差：

ν_l = 2π Δf τ_l

把 s_l 在频率维上截掉最后一个频点，与截掉第一个频点相比，相差 e^{-jν_l}。这就是频率维的旋转不变性。

两个旋转不变性分别对应两个矩阵方程。把信号子空间 Us 按频率-阵元顺序重排后做截断，得到四个子矩阵：

U_θ1 P_θ = U_θ2 U_τ1 P_τ = U_τ2

P_θ 的特征值是 e^{jμ_l}，P_τ 的特征值是 e^{-jν_l}。到这里，角度和时延都变成了两个 L×L 矩阵的特征值问题。但真正麻烦的事情才刚开始：两个矩阵各自做特征分解，特征向量的顺序是任意的，直接按排序对齐极容易配错。

2.3 JADE 的配对手段：相似变换把时延“绑”到角度上

既然 P_θ 和 P_τ 作用在同一个信号子空间 Us 上，它们的特征向量在理想情况下应该张成同一组坐标。只是数值上分别做 eig 分解时，列顺序和归一化方向都不同。我的做法是：先对 P_θ 做特征分解，取特征向量矩阵 V_θ，然后把 P_τ 变换到这个基下：

G = V_θ^{-1} P_τ V_θ

在无噪声情况下，G 是对角阵，对角元素就是配好对的时延旋转量。第 l 个角度特征值对应的时延，就是 G(l,l)。噪声存在时，G 的非对角项不为零但通常远小于对角项，直接取对角元素即可。这一步是 JADE 和“先角度后时延的两步法”最大的区别：配对不是事后做最近邻匹配，而是在同一个子空间基下天然完成的。

参数选择上有几个硬性边界需要记住：阵元数 M 必须大于路径数 L，频点数 Nf 也必须大于 L；否则对应维度的旋转不变方程无法提供足够独立方程。同时，不模糊条件要求空间间距 d≤λ/2、频率间隔 Δf 满足 τ_max < 1/Δf。下表是我后面仿真里会用到的一组参数及其边界依据。

3. 用 MATLAB 搭建仿真链路：多径信道频域快照的生成与参数设定

原理立住之后，先做仿真而不是直接上实采数据，是我一贯的顺序。仿真能让你在已知真值的条件下验证算法每个环节是否正确，尤其是配对那一步。下面这套信号生成代码不依赖任何工具箱，基础 MATLAB 就能跑。

3.1 仿真参数与可复现设定

这里选参数的原则是：让角度和时延都有一定挑战性，但不刻意刁难算法。θ = [-35, 5, 38]° 覆盖了负角度、接近法线和大角度三种情况；τ = [0, 38, 120] ns 中 0 ns 是一条直射径，38 ns 刻意设置得小于常规傅里叶分辨间隔，用来观察超分辨能力。总带宽为 Nf×Δf = 16 MHz，对应的常规时延分辨间隔约 62.5 ns，38 ns 这个值正好在传统分辨率之下，能展示 ESPRIT 类算法的超分辨特性。

角度分辨率主要由阵列孔径决定，8 元半波长间距的 ULA 波束宽度在十几度量级，所以三个角度间隔 30° 以上不算难；真正容易翻车的是角度和时延同时接近时，子空间提取不稳定。后面第 5 章会专门讲这种情况。

3.2 频域多径快照生成代码

先给出完整可运行的信号生成部分。真值、阵列和频响一起定义，最后按 K 个导频符号生成数据矩阵 Y。

% JADE_2D_ESPRIT_demo.m —— 信号生成部分 clear; clc; rng(2024); M = 8; % 阵元数 Nf = 32; % 频域采样点数 df = 0.5e6; % 频率间隔 0.5 MHz d = 0.5; % 阵元间距，以波长为单位 L = 3; % 多径数 K = 200; % 导频符号数（快照数） SNR_dB = 25; % 每径信噪比 % 真实参数 theta_true = [-35, 5, 38]; % 到达角，单位：度 tau_true = [0, 38, 120] * 1e-9; % 时延，单位：秒 % 频点序列 f = (0:Nf-1).' * df; % 阵列导向矢量矩阵 A：M x L A = exp(1j * 2 * pi * d * (0:M-1).' * sind(theta_true)); % 频率维相位移矩阵 E：Nf x L，注意这里指数是负号 E = exp(-1j * 2 * pi * f * tau_true); % 生成 K 个频域快照，每列的复增益 β 随机变化 Y = zeros(M * Nf, K); for k = 1:K beta = exp(1j * 2 * pi * rand(L, 1)); % 每个符号的随机复增益 Xk = A * diag(beta) * E.'; % M x Nf 频域响应 Y(:, k) = Xk(:); % 向量化，按列优先 end % 加入复高斯噪声 noise = (randn(size(Y)) + 1j * randn(size(Y))) / sqrt(2) * 10^(-SNR_dB/20); Y = Y + noise;

这段代码里有一处最容易出错：Xk(:)是按列优先展开的，也就是说向量中前 Nf 个元素属于阵元 1 的全部频点，接着是阵元 2 的全部频点。后面从 Us 重建三维数组时，必须保持同样的排列顺序，否则角度和时延会直接互换。E.'是转置而不是共轭转置，因为频域响应的相位项就是 e^{-j2πfτ}，取复共轭会把符号搞反，导致时延估计出现镜像。

3.3 为什么快照之间必须让复增益变化

很多第一次跑通这段代码的人会试图简化：直接让 beta 固定，只生成一列 Y。结果 SVD 出来的第二奇异值接近零，角度和时延全都估计不出来。原因在于，如果每个快照的复增益完全相同，数据矩阵 Y 的每一列都是同一个联合导向矢量的线性组合，秩最多为 1，信号子空间根本张不起来。

我一般会让 beta 的相位在 0 到 2π 之间随机变化，相当于模拟 OFDM 符号中不同数据调制带来的相位变化。实际信道探测里，发射端连续发送不同的伪随机序列或 OFDM 符号，接收端自然就能得到这种变化。K 取 200 是经验值，保证信号子空间在 25 dB 噪声下有足够统计稳定性；K 太小时奇异值展开不充分，配对容易受噪声扰动。

4. 核心实现：二维 ESPRIT 与 JADE 自动配对的可运行代码

信号生成完毕，接下来是核心的估计部分。整个过程分三步：SVD 提取信号子空间、用截断矩阵构建两个旋转算子、用相似变换完成配对并还原角度和时延。

4.1 从数据矩阵到信号子空间

先对 Y 做 SVD，取前 L 个主左奇异向量作为信号子空间 Us。这里的 L 可以先按真实路径数给，第 6 章会讲如何用信息论准则自动定阶。

% 第 1 步：SVD 提取信号子空间 [U, S, ~] = svd(Y, 'econ'); Us = U(:, 1:L); % MNf x L 信号子空间 sv = diag(S); % 观察奇异值分布，确认秩 figure; plot(10*log10(sv), 'o-'); xlabel('奇异值序号'); ylabel('dB'); title('奇异值分布，前 L 个应明显高于噪声平台');

这里有个判断信号子空间提取是否正常的技巧：画出奇异值分布，如果前 L 个奇异值与后面的噪声平台有 20 dB 以上的落差，说明数据秩充分、信号模型匹配。如果第二奇异值掉得特别快，通常就是快照间 β 没有变化；如果前 L 个之后还有接近的平台，可能是实际路径数大于设定的 L。

4.2 两个平移不变方程的最小二乘求解

构建旋转算子的关键在于把 Us 重排成三维数组，并保证维度顺序与向量化时一致。正如 3.2 中说的，向量里前 Nf 个元素是阵元 1 的频点，所以重排时必须把 Nf 放在第一维，M 放在第二维。

% 第 2 步：重排信号子空间，构造平移不变方程 Umat = reshape(Us, [Nf, M, L]); % 第一维：频点，第二维：阵元 % 角度维：相邻阵元之间的旋转不变性 Uth1 = reshape(Umat(:, 1:M-1, :), [], L); % 去掉第 M 个阵元 Uth2 = reshape(Umat(:, 2:M, :), [], L); % 去掉第 1 个阵元 Pth = pinv(Uth1) * Uth2; % L x L 空间旋转算子 % 时延维：相邻频点之间的旋转不变性 Uta1 = reshape(Umat(1:Nf-1, :, :), [], L); % 去掉最后频点 Uta2 = reshape(Umat(2:Nf, :, :), [], L); % 去掉首频点 Pta = pinv(Uta1) * Uta2; % L x L 频率旋转算子

pinv(Uth1) * Uth2是标准最小二乘解。当噪声较大时，这个解略有偏置，因为等号两侧的 Uth1 和 Uth2 都含噪声。更稳的做法是用总体最小二乘，但先跑通最小二乘版本、确认配对逻辑正确，再升级也不迟。Umat 的维度顺序必须有注释写清楚，这是最容易翻车的地方，后面第 5 章会专门讲。

4.3 相似变换配对与角度/时延还原

得到 Pth 和 Pta 之后，不能直接分别做 eig 再按大小排序。因为两个矩阵的特征向量虽然理论上张成同一组方向，但数值分解时列顺序不同，直接对齐必然配错。正确做法是取 Pth 的特征向量矩阵作为公共基，把 Pta 变换过去。

% 第 3 步：用相似变换完成联合配对 [Vth, Dth] = eig(Pth); % Vth：角度特征向量，Dth：角度特征值 G = Vth \ (Pta * Vth); % Pta 在 Vth 基下的表示 % 提取旋转量 mu = angle(diag(Dth)); % 空间旋转相位 e^{jμ} nu = angle(diag(G)); % 频率旋转相位 e^{-jν} % 还原角度和时延 theta_est = asind(mu / (2 * pi * d)); tau_est = -nu / (2 * pi * df); % 负号来自 e^{-j2πfτ} % 按角度排序显示 [theta_est, idx] = sort(theta_est); tau_est = tau_est(idx); disp('估计角度 (deg) vs 真值:'); disp([theta_est(:), theta_true(:)]); disp('估计时延 (ns) vs 真值:'); disp([tau_est(:)*1e9, tau_true(:)*1e9]);

angle(diag(Dth))取的是特征值复角的相位主值，范围在 -π 到 π 之间。空间中 μ = 2πd sinθ，只要 d≤0.5λ，这个范围内不会出现角度模糊；频率维 ν = 2πΔfτ，只要 τ<1/Δf，同样不会出现时延折叠。tau_est前面的负号是沿自频域相位项 e^{-j2πfτ}，漏掉负号会让时延估计变成负值或镜像值。

配对逻辑需要多说一句：G 是对角占优矩阵，第 l 个对角元素恰好就是与第 l 个角度特征值相配的时延旋转量。这一特性来自 Pth 和 Pta 共享信号子空间这一事实。如果噪声很大，G 的非对角项变大，diag(G) 不再严格对应，这是第 5 章要单独处理的。

5. 避坑与参数调优：我在 JADE MATLAB 代码上踩过的 5 个坑

这套代码我在不同项目里跑过很多遍，每次翻车几乎都集中在少数几个点上。有些坑光看报错根本看不出来，输出数值看着正常但实际是错的。下面这 5 条是最常见的。

5.1 把 JADE 当成 ICA 的联合对角化，方向全错

现象：搜索 JADE 相关程序时，下载下来的代码输出的是“源信号”矩阵，完全没有角度和时间的概念，运行后得到的东西和阵列测向毫无关系。

原因：信号处理领域里有两个 JADE。一个是本篇讲的 Joint Angle and Delay Estimation，另一个是独立成分分析里的 Joint Approximate Diagonalization of Eigenmatrices。后者在 MATLAB 生态里更常见，很多名为 jade.m 的程序都是 ICA 工具。

解决：下载代码时先看函数输入输出。ICA 的 JADE 输入是观测信号矩阵，输出是解混矩阵和源信号；本标题的 JADE 输入是阵列频域快照，输出是角度与时延估计。文件名建议直接写JADE_2D_ESPRIT_demo.m，避免混淆。另外，ICA 版本的 JADE 一般不会出现阵列导向矢量和频域相位项，只要代码里有这两个概念，才是正确的方向。

5.2 快照间复增益不变，子空间秩退化

现象：SVD 之后第二个奇异值比第一个小几个数量级，估计出的角度和时延要么只有一个值，要么完全随机。

原因：3.3 里已经说过，如果每个导频符号的 β 都相同，Y 的每一列都是同一组导向矢量的固定线性组合，数据矩阵的秩退化到接近 1。常见诱因是仿真代码把 beta 定义在循环外面，或者实测数据里导频符号完全重复。

解决：确认每列生成时 beta = exp(1j2pi*rand(L,1)) 在循环内部；实测场景中，如果导频符号之间相关性太强，可以人为在所有导频符号上乘以一个已知的随机相位序列，只要接收端知道这个序列，解调后 β 的变化就保留了。

5.3 reshape 按列优先摆错维度，角度和时延互换

现象：估计出的数值范围倒是合理，但看起来像是 θ 和 τ 互换了，比如 θ 估计出几个零点几的值，τ 却变成几十度。

原因：Xk(:)按列优先展开，向量中前 Nf 个元素属于阵元 1 的全部频点。如果重排时写成reshape(Us, [M, Nf, L])，第一维变成阵元、第二维变成频点，构造角度和时延平移方程时截断的维度就反了。

解决：一律按reshape(Us, [Nf, M, L])重排，并在代码里用注释锁定这个约定。跑通后可以用一个最简模型自检：设 θ_true = [0, 10]°，τ_true = [0, 50] ns，如果估计出来的 θ 和 τ 数量级正确，说明维度没问题。这个自检 5 分钟就能做完，比反复调试代码快得多。

5.4 时延超过 1/Δf，估计值折叠到错误区间

现象：真实时延是 2.3 μs，估计结果却落在 0.3 μs 附近，且改变仿真种子后数值在某个区间内跳变。

原因：频率旋转量 ν = 2πΔfτ 以 2π 为周期，angle() 取主值后只能表示 [-π, π]，对应 τ 的可测范围是 [-1/(2Δf), 1/(2Δf)]。超过这个范围，真实时延就被折叠回主值区间。

解决：先根据场景先验压低 Δf，让 1/Δf 大于最大可能时延。实际信道探测中，如果时延扩展很大，可以先用相关法做一个粗时延对齐，把目标时延搬移到不模糊区间内，再做 JADE 精估计。这个“粗对齐 + 精估计”的组合比单纯增大 Nf 更经济。

5.5 噪声上来之后配对乱跳，先看 G 的非对角项

现象：SNR 低于 15 dB 时，固定种子偶尔正确，换一个种子就出现角度对但时延错，或者整体配对乱序。

原因：噪声增大后 Pth 和 Pta 自身的特征向量扰动变大，G 不再是对角占优矩阵，diag(G) 无法可靠代表匹配结果。这时最直观的现象就是 G 的非对角项与对角项在同一数量级。

解决：先看abs(G)的矩阵分布，确认是否对角占优。如果不占优，先把 Pth 和 Pta 的求解从最小二乘换成总体最小二乘，通常能改善 3-5 dB 的可用信噪比；再把 K 从 200 提高到 400 或 800，用更多快照压低子空间估计方差。如果这两个手段都不够，就要怀疑路径数 L 定得不对，回 4.1 检查奇异值分布。

6. 进阶用法：从仿真走向实测前的三个习惯

仿真跑通之后，真正把它用在实际数据上，我建议先养成三个习惯。第一个习惯是自动定阶：不要永远手填 L。实测中多径数未知，我一般会先看奇异值，再用信息论准则如 MDL 自动估计路径数。做法是把奇异值按降序排列，计算噪声特征值的算术平均值与几何平均值之差，路径数通常取在两者最接近的位置附近。这一步不用写得很严格，核心是让程序对真实场景的路径数变化有一定的容忍度。

第二个习惯是适配 OFDM 实测数据。实际的宽带通信系统里，接收端经过信道估计后能得到每个 OFDM 符号在每个子载波上的频域响应，这正好就是第 3 章数据矩阵 Y 的列。每个 OFDM 符号的数据调制天然提供了 β 的随机变化，所以只需要把 Y 替换成 OFDM 信道估计结果，把 f 换成对应的子载波频率序列，算法主体不用改动。需要注意两点：一是剔除掉导频空洞或零子载波对应的频点，不要让频率间隔变成非均匀的；二是如果信道随时间慢变，K 个符号的采集时间不能超过信道相干时间，否则 β 的变化中混入了信道本身的衰落，子空间会被污染。

第三个习惯是验证时不要只看一次估计值。我会固定真值，在 SNR = 0 到 30 dB 之间循环仿真，计算角度和时延估计的 RMSE，画成随 SNR 变化的曲线。这条曲线会有一个明显的阈值效应：高于某个 SNR 时误差迅速下降，低于阈值时配对开始乱跳。这个阈值就是这套参数下算法的可用下限。如果你要部署到实际系统里，先用这条路找清楚阈值，比在实采数据上反复试参数要高效得多。

最后提一个我自己的习惯：任何新场景里，我都会跑一次 θ_true = [0, 10°]、τ_true = [0, 50 ns] 的最小自检，确认维度布局没有因为数据排列方式变化而出错。排摸信号处理的问题，先确认模型与代码布局一致，再去怀疑噪声和配对，这个顺序能省掉一半的排查时间。希望这篇笔记能帮你在自己的阵列上也跑出稳定的二维参数估计。

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