简介:基于正则化约束总体最小二乘(RCTLS)的单站DOA-TDOA无源定位算法复现文档,面向无线通信、雷达、导航定位领域的科研人员与工程师,旨在解决单站联合利用到达角度(DOA)与到达时差(TDOA)进行定位时面临的方程病态和精度不足问题。资源为单个docx文档,压缩包共1个文件,大小仅53KB,但内容涵盖完整可运行Python代码及逐步解释,包括RCTLS_DOATDOA类初始化、观测方程线性化、代价函数构建、牛顿迭代求解器、最优正则化参数自动选择以及几何精度因子(GDOP)计算等模块。文档还深入讨论了系数矩阵病态问题的理论分析与改进实现,给出低高度目标定位、误差建模、鲁棒求解器优化等扩展方向,并通过仿真实验验证RCTLS相比传统CTLS在定位精度和鲁棒性上的优越性,同时分析了目标与辐射源位置对定位精度的影响。已有73人学习浏览,适合需要理论推导与工程实现并重,或希望将改进算法用于单站无源定位系统及外辐射源布局优化的开发者参考。
1. 单站 DOA-TDOA 无源定位为什么需要 RCTLS
做单站无源定位的同行应该都有这种体会:测向方程和时差方程各自都能列出来,但联立之后解出来的目标点总是不对劲,尤其是目标落在反射路径和直达路径夹角很小的区域时,误差能到几百米。问题不在测向仪或时差检测器上,而在求解方式上——你用的最小二乘默认系数矩阵是精确的,可这个场景里 DOA 测量误差已经顺着三角函数爬进了系数矩阵。RCTLS(正则化约束总体最小二乘)就是为这种情况设计的:它同时承认系数矩阵和观测向量都带噪声,再用正则化项压制病态方程的解发散。这篇文章按论文复现的思路,从公式推导讲到 MATLAB 实现,再到参数调整和验证方法,适合正在复现定位类论文、或者被测向测时差数据折磨的工程师照着落地。
2. RCTLS 的理论拆解:为什么 LS 会在 DOA-TDOA 上失效
2.1 把 DOA 和 TDOA 写成 Ax=b:方向线与路径差的两个线性化
单站 DOA-TDOA 定位的模型假设是:一个观测站接收目标信号,除了直达径之外还有来自反射体的多径信号。观测站测到两个量——目标信号的到达方位角,以及直达径与反射径之间的到达时间差。观测站位置、反射点位置已知,目标坐标是待求量。
先看 DOA。设观测站坐标为 $x_s=[x_{s},y_{s}]^T$,目标位置为 $p=[p_x,p_y]^T$,方位角 $\alpha$ 定义为从 x 轴逆时针转到目标方向的角度。目标在方向线上,等价于目标相对观测站的矢量与方向单位矢量 $u=[\cos\alpha,\sin\alpha]^T$ 平行,也就是法向分量等于零。取 $u_\perp=[-\sin\alpha,\cos\alpha]^T$,DOA 方程就是:
$u_\perp^T(p-x_s)=0$
展开后得到一个关于 p 的线性方程:
$-\sin\alpha\cdot p_x+\cos\alpha\cdot p_y=-\sin\alpha\cdot x_s+\cos\alpha\cdot y_s$
这一行里 $\alpha$ 是测量值,所以系数 $-\sin\alpha$ 和 $\cos\alpha$ 都带测量噪声,右边也是。系数矩阵带噪声这件事,就是 TLS 和 RCTLS 登场的理由。
TDOA 部分要绕一个弯。第 i 条反射路径对应的反射点记为 $x_{r,i}$,测得的直达径与反射径时差为 $\tau_i$,统一转成距离差 $r_i=c\tau_i$。直达径长度等于观测站到目标的距离 $d=|p-x_s|$,于是反射径长度可以写成 $d-r_i$(符号约定后面避坑章会专门讲)。由几何关系,反射径长度又等于 $|p-x_{r,i}|$,而 $p=x_s+d\cdot u$,代入并消去 $d^2$ 项后得到:
$d=\frac{r_i^2-|x_s-x_{r,i}|^2}{2\left(u^T(x_s-x_{r,i})+r_i\right)}$
这个式子的含义是:每条 TDOA 路径都能独立给出一个目标到观测站的距离估计 $d_i$,$d_i$ 的噪声同时来自 $\tau_i$ 的测量误差和由 $\alpha$ 引入的 $u$ 的误差。把它写成纵向约束:
$u^T(p-x_s)=d_i$
现在把 DOA 的横向约束和每条 TDOA 的纵向约束叠起来,就得到一个超定线性系统 $Jp=y$。例如一个 DOA 加两条反射路径时:
$J=\begin{bmatrix}-sin\alpha & cos\alpha \ cos\alpha & sin\alpha \ cos\alpha & sin\alpha\end{bmatrix}, \quad y=\begin{bmatrix}u_\perp^T x_s \ u^T x_s+d_1 \ u^T x_s+d_2\end{bmatrix}$
这个形式的好处是,所有测量误差都集中体现在 $J$ 和 $y$ 上,后面无论用 LS、TLS 还是 RCTLS,输入都是同一对矩阵。
2.2 TLS 和 LS 的分水岭:系数矩阵噪声不能假装不存在
最小二乘的隐含假设是系数矩阵精确、只有观测向量有噪声。求解 $Jp=y$ 时,LS 给出的是让 $|Jp-y|2^2$ 最小的解,闭式解为 $p{LS}=(J^TJ)^{-1}J^Ty$。
问题是这个模型里 $J$ 的每一行都是拿带噪的 $\alpha$ 算出来的三角函数。DOA 误差 0.5 度时,$\sin\alpha$ 的相对误差在 0.0087 量级,看起来不大,但当 $J$ 接近病态、条件数达到 $10^4$ 以上时,这个误差会被放大到完全不可接受的程度。LS 的解在此情况下是有偏的,而且偏差方向不固定,跟目标相对观测站的方位相关。
总体最小二乘把噪声假设往前推了一步:系数矩阵和观测向量都有噪声。它的目标是找一个尽可能小的扰动 $[\Delta J,\Delta y]$,使得扰动后的方程组 $(J+\Delta J)p=y+\Delta y$ 严格相容。这个目标可以等价写成瑞利商形式:
$\min_p \frac{|Jp-y|^2}{1+|p|^2}$
分母里的 $|p|^2$ 来自对 $\Delta J$ 的 Frobenius 范数最小化,这不是人为加进去的,而是 TLS 问题的数学结构。求这个最小值不需要迭代优化,对增广矩阵 $[J,y]$ 做 SVD,取最小奇异值对应的右奇异向量 $v$,令 $p_{TLS}=v(1:n)/v(n+1)$ 即可。注意分母是 $v$ 的最后一个分量,如果它接近零,说明问题退化,TLS 解会非常不稳定,这一刻就开始需要正则化了。
2.3 正则化项如何压制病态:RCTLS 的目标函数与梯度
TLS 解决了系数矩阵噪声的问题,但在单站 DOA-TDOA 场景里还不够。反射点、观测站、目标三者之间的几何关系一旦接近退化,$J$ 的条件数会急剧膨胀。比如目标恰好落在观测站与反射点的延长线附近时,TDOA 方程和 DOA 方程几乎线性相关,$J$ 的最小奇异值趋近于零,TLS 解虽然理论上无偏,但方差大得没有实用价值。
正则化的思路是给解的范数加一个惩罚项,压制病态带来的解漂移。RCTLS 的目标函数写为:
$F(p)=\frac{|Jp-y|^2}{1+|p|^2}+\lambda|p|^2$
这里 $\lambda$ 是正则化参数。$\lambda=0$ 时退化为 TLS;把分母去掉则退化为岭回归,但那样就丢失了系数噪声校正能力。这个目标函数第一项修正系数矩阵噪声,第二项控制解的范数,两者缺一不可。从约束的角度看,正则化项等价于给解空间加了一个先验约束——这也就是"约束"两个字在 RCTLS 里的含义。
求解这个函数我采用解析梯度加回溯线搜索。设 $r=Jp-y$,$N=r^Tr$,$D=1+p^Tp$,梯度为:
$\nabla F=\frac{2J^TrD-2Np}{D^2}+2\lambda p$
这个函数不是凸的,所以初值选择很关键,第 4 章会专门讲初值策略。$\lambda$ 较大时,正则项给目标函数增加了凸性,局部极小点会变少,但 $\lambda$ 太大会把解拉向原点,变成有偏估计,这个矛盾是所有正则化方法都要面对的。
3. 用 MATLAB 复现 RCTLS 单站定位:从场景生成到三种方法对比
3.1 场景与测量生成:先跑通单次定位再谈统计
复现定位算法第一步不是写求解器,而是把场景参数定下来。这里用一个典型的多径单站场景:观测站在坐标原点,两个反射点分别模拟建筑物墙体和地面反射,目标放在大约 3 公里外。这个尺度下,DOA 误差 0.5 度对应的横向误差约为 26 米,TDOA 误差 5 纳秒对应的距离误差约为 1.5 米,量级差异明显,正好能看出加权处理的价值。
% 场景与测量生成:单站多径 DOA-TDOA xs = [0; 0]; % 观测站位置 (m) xr1 = [800; 300]; % 反射点 1(墙面)(m) xr2 = [-400; 900]; % 反射点 2(地面/建筑)(m) p_true = [3000; 2200]; % 目标真实位置 (m) c = 3e8; % 光速 (m/s) % 真实方向角与真实时差 u_true = (p_true - xs) / norm(p_true - xs); alpha_true = atan2(u_true(2), u_true(1)); tau1_true = (norm(p_true - xs) - norm(p_true - xr1)) / c; tau2_true = (norm(p_true - xs) - norm(p_true - xr2)) / c; % 添加测量噪声 sigma_alpha = 0.5 * pi / 180; % 测向标准差 0.5 度 sigma_tau = 5e-9; % 时差标准差 5 ns alpha = alpha_true + sigma_alpha * randn; tau1 = tau1_true + sigma_tau * randn; tau2 = tau2_true + sigma_tau * randn; % 由带噪测量构造线性系统 J * p = y u = [cos(alpha); sin(alpha)]; r1 = c * tau1; r2 = c * tau2; s1 = u' * (xs - xr1); s2 = u' * (xs - xr2); L1 = (xs - xr1)' * (xs - xr1); L2 = (xs - xr2)' * (xs - xr2); % 每条路径独立估计目标到观测站的距离 d1 = (r1^2 - L1) / (2 * (s1 + r1)); d2 = (r2^2 - L2) / (2 * (s2 + r2)); J = [-sin(alpha), cos(alpha); % DOA 横向约束 cos(alpha), sin(alpha); % TDOA 路径 1 纵向约束 cos(alpha), sin(alpha)]; % TDOA 路径 2 纵向约束 y = [-sin(alpha)*xs(1) + cos(alpha)*xs(2); cos(alpha)*xs(1) + sin(alpha)*xs(2) + d1; cos(alpha)*xs(1) + sin(alpha)*xs(2) + d2];这段代码里最需要注意的是 $d_i$ 的计算。$\tau_i$ 的符号定义是直达径到达时间减去反射径到达时间,所以 $r_i$ 是正的,分母 $s_i+r_i$ 的正负取决于反射点相对观测站和目标的几何关系。单次定位如果跑到镜像位置,先检查这个符号,不必急着怀疑求解器。
3.2 RCTLS 求解器:解析梯度与回溯线搜索的实现
RCTLS 求解器是整个复现的核心。我不用现成的 fminunc,而是手写梯度下降加回溯线搜索,这样每一步都在掌控之中,也方便在论文里解释算法细节。
function x = rctls_solve(J, y, lambda, opts) % RCTLS 求解器:迭代求解正则化总体最小二乘 % 目标函数: F(x) = ||Jx-y||^2 / (1+||x||^2) + lambda * ||x||^2 % 输入: % J : m x n 带噪系数矩阵 % y : m x 1 带噪观测向量 % lambda : 正则化参数 % opts : 可选结构体,max_iter / tol / verbose % 输出: % x : 定位结果 if nargin < 4, opts = struct(); end max_iter = 200; tol = 1e-10; if isfield(opts, 'max_iter'), max_iter = opts.max_iter; end if isfield(opts, 'tol'), tol = opts.tol; end % 初值使用岭回归解,保证在病态下也有合理起点 [m, n] = size(J); x = (J'*J + lambda*eye(n)) \ (J'*y); for iter = 1:max_iter r = J*x - y; N = r'*r; D = 1 + x'*x; F_old = N/D + lambda * (x'*x); % 解析梯度 grad = (2*J'*r*D - 2*N*x) / D^2 + 2*lambda*x; % 回溯线搜索:保证充分下降 alpha = 1.0; while alpha > 1e-12 x_new = x - alpha * grad; r_new = J*x_new - y; F_new = (r_new'*r_new)/(1 + x_new'*x_new) + lambda*(x_new'*x_new); if F_new <= F_old - 1e-4 * alpha * (grad'*grad) break; end alpha = alpha * 0.5; end if norm(x_new - x, inf) < tol x = x_new; break; end x = x_new; end end这个实现里有几个参数值得说。充分下降准则里的 $10^{-4}$ 是优化教材里的标准值,太小会让线搜索退化成全步长,太大则步长衰减过快。回溯线搜索的初始步长设为 1,如果目标函数不够光滑,单次梯度下降可能跳出合理区域,这时线搜索会自动把步长减半。迭代终止条件用的是无穷范数小于 $10^{-10}$,对坐标量级在千米的定位问题来说,这个阈值已经远低于双精度浮点的有效位数了。
初值选择用岭回归解而不是纯 LS 解,是因为当 $J$ 条件数很差时,$(J^TJ)^{-1}$ 本身就不稳定,加一个小的 $\lambda I$ 至少保证矩阵可逆,且解在合理范围内。
3.3 与 LS / TLS 对比:Monte Carlo 循环与误差统计
单次定位说明不了问题,RCTLS 的论文复现一定需要 Monte Carlo 统计。TLS 的 SVD 解法实现很简单,可以直接作为对比基准:
function x = tls_solve(J, y) % 总体最小二乘:对增广矩阵做 SVD,取最小奇异值右奇异向量 [~, ~, V] = svd([J, y]); v = V(:, end); if abs(v(end)) < 1e-12 x = J \ y; % 退化时退化为 LS else x = v(1:end-1) / v(end); end end主程序里对三种方法跑同样的随机噪声,统计 RMSE:
% Monte Carlo 对比 LS / TLS / RCTLS N = 500; err_ls = zeros(N,1); err_tls = zeros(N,1); err_rctls = zeros(N,1); lambda = 1e-6 * norm(J, 'fro')^2; % 正则化参数经验量级,第 4 章详述 for k = 1:N alpha = alpha_true + sigma_alpha * randn; tau1 = tau1_true + sigma_tau * randn; tau2 = tau2_true + sigma_tau * randn; u = [cos(alpha); sin(alpha)]; r1 = c*tau1; r2 = c*tau2; s1 = u'*(xs-xr1); s2 = u'*(xs-xr2); d1 = (r1^2-L1)/(2*(s1+r1)); d2 = (r2^2-L2)/(2*(s2+r2)); Jk = [-sin(alpha), cos(alpha); cos(alpha), sin(alpha); cos(alpha), sin(alpha)]; yk = [-sin(alpha)*xs(1)+cos(alpha)*xs(2); cos(alpha)*xs(1)+sin(alpha)*xs(2)+d1; cos(alpha)*xs(1)+sin(alpha)*xs(2)+d2]; p_ls = Jk \ yk; p_tls = tls_solve(Jk, yk); p_rctls = rctls_solve(Jk, yk, lambda); err_ls(k) = norm(p_ls - p_true); err_tls(k) = norm(p_tls - p_true); err_rctls(k) = norm(p_rctls - p_true); end fprintf('LS RMSE: %.2f m\n', sqrt(mean(err_ls.^2))); fprintf('TLS RMSE: %.2f m\n', sqrt(mean(err_tls.^2))); fprintf('RCTLS RMSE: %.2f m\n', sqrt(mean(err_rctls.^2)));把这段代码里的场景参数、噪声参数换成你自己论文里的配置,就能复现大部分单站定位实验。代码讲解要连数据流一起看:$J$ 和 $y$ 的每一行都来自带噪测量,所以 LS 的假设从一开始就不成立,TLS 和 RCTLS 的优势在 Monte Carlo 统计里会体现为 RMSE 的明显下降。
4. RCTLS 的 4 个必调参数:正则化系数、加权矩阵与收敛设置
4.1 正则化参数 λ:L 曲线选点与经验量级
$\lambda$ 是 RCTLS 里最敏感的参数。实际调试时,我会先固定一个场景,扫描 $\lambda$ 从 $10^{-8}$ 到 $10^2$ 的取值,画出 $|Jp-y|$ 对 $|p|$ 的曲线——也就是 L 曲线。$\lambda$ 偏小时解范数大,残差小但方差大;$\lambda$ 偏大时解范数被压住,残差增大,出现明显偏差。L 曲线的拐角处曲率最大,那一点对应的 $\lambda$ 就是折中值。
% L 曲线选点:扫描 lambda,画 residual-norm 对 solution-norm lambdas = logspace(-8, 2, 30); resid_norm = zeros(size(lambdas)); sol_norm = zeros(size(lambdas)); for i = 1:numel(lambdas) x = rctls_solve(J, y, lambdas(i)); resid_norm(i) = norm(J*x - y); sol_norm(i) = norm(x); end loglog(sol_norm, resid_norm, 'o-');如果不想每次都画图,有一个省事的经验公式:$\lambda_0 = 10^{-6}\cdot|J|_F^2$,其中 $|J|_F$ 是 Frobenius 范数。这个量级在目标距离数千米、$J$ 元素为 1 量级的场景下表现稳定。注意这个公式只是起点,正式做实验前还是要用 L 曲线确认一遍,因为几何退化程度不同,最优 $\lambda$ 会差几个数量级。
4.2 加权矩阵:让 DOA 噪声和 TDOA 噪声在同一个量纲下比较
我在第 3 章的代码里没有加权,但实际跑出来的结果往往有一条规律:TDOA 测得很准(纳秒级),DOA 测得不差(零点几度),但定位误差仍然很大。原因是两者对最终解的影响量级不匹配。0.5 度的 DOA 误差在 3 公里距离上折算成横向距离误差约 26 米,而 5 纳秒的 TDOA 误差只折算成 1.5 米。如果直接把两行方程叠在一起,TDOA 这条信息几乎被 DOA 噪声主导,相当于只用了方向没用到时差。
解决办法是加权。给每个方程行的残差除以该行的等效噪声标准差,等价于在目标函数里引入对角权重矩阵 $W$:
$F(p)=\frac{|W(Jp-y)|^2}{1+|p|^2}+\lambda|p|^2$
实现时不必改求解器,只要对 $J$ 和 $y$ 左乘 $W$ 再调用 rctls_solve 就行。权重按经验设:DOA 行的等效距离噪声 $\sigma_{eq}=d\cdot\sigma_\alpha$,TDOA 行用 $c\cdot\sigma_\tau$。$d$ 可以用一次 LS 解先估出来,代入权重后再跑 RCTLS。这个两遍处理的做法在论文里也常见,不算偷懒。
4.3 迭代、收敛与初值:RCTLS 的数值稳定性设置
RCTLS 目标函数非凸,初值决定了最后收敛到哪个局部极小点。我一般把关卡设在初值上:第一优先级是岭回归初值,它比 LS 初值在病态时更稳;第二优先级是四方向重试,即从 $x_0$、$x_0+[\Delta,0]^T$、$x_0+[0,\Delta]^T$、$x_0-[\Delta,\Delta]^T$ 四个点分别迭代,取最终目标函数值最小的解。这里的 $\Delta$ 取场景尺度的千分之一,比如场景 3 公里就取 3 米。这个技巧在目标接近几何退化区域时特别管用。
回溯线搜索的参数我固定用充分下降系数 $10^{-4}$、步长衰减因子 0.5。最大迭代次数 200 次足够让梯度降到 $10^{-12}$ 量级,迭代终止阈值 $10^{-10}$ 对坐标解来说已经是过度收敛了,主要作用是防止在极小点附近做无用功。如果你的场景坐标是经纬度(量级 $10^6$),记得把这两个阈值按比例放大,否则收敛判据可能被浮点误差卡住。
4.4 一张参数表:从 λ 到 Monte Carlo 次数
| 参数 | 推荐值 | 说明 |
|---|---|---|
| $\sigma_\alpha$ | 0.1°~2° | 典型阵列测向精度,小于 0.1° 时噪声基本由多径衰落主导 |
| $\sigma_\tau$ | 1~50 ns | 相关峰值时差估计,受带宽和信噪比影响 |
| $\lambda$ | $10^{-6}||J||_F^2$ | L 曲线拐角,最大不超过 $J$ 最大奇异值的 $10^{-2}$ 倍 |
| 加权矩阵 $W$ | 每行除以其等效噪声标准差 | DOA 行用 $d\sigma_\alpha$,TDOA 行用 $c\sigma_\tau$ |
| max_iter | 200 | 梯度降到 $10^{-12}$ 所需迭代次数的两倍安全量 |
| tol | $10^{-10}$ | 坐标量级为 $10^3$ 时的经验值,按场景坐标量级缩放 |
| Monte Carlo 次数 | 200~500 | 少于 100 次时 RMSE 抖动明显,500 次以上曲线稳定 |
这张表不是拍脑袋,是拿第 3 章代码跑了不同信噪比组合后得到的行为规律。$\lambda$ 和 $W$ 是唯一需要逐场景调整的两项,其余参数一套配置可以用到整个实验序列里。
5. RCTLS 复现翻车记录:5 个让结果崩掉的问题与排查
5.1 现象一:DOA 误差在小偏置几何下被放大成百米级偏差
有次我把反射点放在离观测站只有 200 米的位置,目标在 3 公里外,RCTLS 单次定位误差直接飙到 400 米,而换一个反射点布局就恢复正常。查了很久才发现是 $d_i$ 计算公式里的分母 $s_i+r_i$ 出了问题。$s_i=u^T(x_s-x_{r,i})$ 是观测站到反射点矢量在目标方向上的投影,反射点离观测站很近时,这个投影和 $r_i$ 几乎大小相等、符号相反,分母趋近于零,DOA 测量误差在这个除法中被放大了一个数量级以上。
解决方式很简单:在计算 $d_i$ 之前检查 $|s_i+r_i|$,如果小于某个阈值(比如场景尺度的 $10^{-3}$ 倍),就丢弃这条 TDOA 路径,只用 DOA 和另一条 TDOA 定位。这不是算法退步,而是明确告诉模型:这条路径的几何条件已经差到测量值无法提供有效信息的程度。
5.2 现象二:TDOA 符号反了,目标直接跑到镜像位置
$\tau$ 的定义在不同论文里不统一。有的定义直达径到达时间减去反射径到达时间,有的反过来。我在复现早期就吃过亏:第 2.1 节的推导默认的是正值(直达径先到),用的是 $r_i=c\tau_i$,反射径长度是 $d-r_i$。如果代码里 $\tau$ 符号相反,$r_i$ 变成负值,$d_i$ 计算出的距离会错误地偏大,定位结果跑到相对反射点的镜像位置。
排查方法很简单:先用真实无噪的 $\alpha$ 和 $\tau$ 跑一遍单次定位,如果误差不是零,先看是不是符号问题,而不是求解器问题。这个自检步骤应该写进主程序开头,堪称后悔药。
5.3 现象三:λ 取大了,RMSE 反而比 LS 还差
这是最让人沮丧的情况:RCTLS 是改进算法,结果不如 LS,那论文怎么写得下去。有一次我把 $\lambda$ 固定成 1,目标距离 3 公里,解被正则项拉向原点,定位结果系统性偏近,Monte Carlo 的 RMSE 比 LS 还大 30%。看 L 曲线才明白,这个场景的最优 $\lambda$ 在 $10^{-5}$ 附近,1 已经属于强正则化区间。
解决思路是不要把 $\lambda$ 当成一个固定常数去调,而是每换一组噪声参数就重新画一次 L 曲线。RCTLS 里 $\lambda$ 的量级由 $|J|$ 决定,噪声水平变了,$\lambda$ 的合理区间也跟着变。给 $\lambda$ 设置一个与 $|J|_F^2$ 成比例的基准值,再在这个基准值附近做小范围扫描,比直接手调稳定得多。
5.4 现象四:量纲不匹配让奇异值分解直接失效
DOA 方程的量级是 1(三角函数值),TDOA 方程转成距离后量级是 $10^3$ 甚至 $10^4$。把它们直接叠进同一个 $J$ 后,SVD 分解和梯度计算会本能地优先拟合那些数值大的行,小量级的 DOA 约束几乎被忽略。TLS 的最小奇异值向量也被带偏,RCTLS 的正则化项实际惩罚的只剩大数坐标。
解决方式是尺度归一化。把距离量纲的行除以一个特征长度 $D_0$,比如观测站到场景中心的粗略距离,解出 $p$ 后再乘回 $D_0$。这样每个方程行的量级都归一到 1 附近,SVD 分解和梯度下降才真正对所有方程一视同仁。第 3 章代码里 $y$ 的量级混着 $10^3$ 和 1,跑通可以,做精细实验前必须先做这一步。
5.5 现象五:结果依赖初值,换一个场景就收敛到局部极小
RCTLS 目标函数非凸,局部极小点的存在和 $\lambda$ 大小相关。$\lambda$ 很小的时候,目标函数接近纯 TLS 的瑞利商,非凸性明显;$\lambda$ 很大时又变成强凸,但解偏差也大。我在某个场景下用岭回归初值能收敛到好结果,换一个反射点布局后同一个初值策略就失效了,定位误差大得离谱。
解决方式是 4.3 节提过的四方向重试。四个初值分别迭代,取目标函数值最小的解。这个策略把非凸优化的风险降低了一个量级,代价是计算量乘四倍。单站定位场景的方程规模只有 3×2,四倍计算量微乎其微,但收敛稳定性收益很大。
6. 性能验证的最后一公里:CRLB、RMSE 曲线与误差椭圆
6.1 用 CRLB 判断算法距离理论下界有多远
复现论文时如果只给三条 RMSE 曲线,审稿人和你自己都会觉得少了点什么。把 CRLB(克拉美-罗界)画上去,才能说明 RCTLS 的性能余量。CRLB 计算需要观测模型关于目标位置的雅可比矩阵,数值差分就足够:
% 数值法计算 CRLB:观测函数 h(p) = [alpha; tau1; tau2] h = @(p) [atan2(p(2)-xs(2), p(1)-xs(1)); (norm(p-xs)-norm(p-xr1))/c; (norm(p-xs)-norm(p-xr2))/c]; H = zeros(3, 2); dp = 1.0; % 差分步长,按坐标量级取 for j = 1:2 e = zeros(2,1); e(j) = dp; H(:, j) = (h(p_true+e) - h(p_true-e)) / (2*dp); end W = diag([1/sigma_alpha^2, 1/sigma_tau^2, 1/sigma_tau^2]); FIM = H' * W * H; crlb = trace(inv(FIM));把 CRLB 开根号后和 RMSE 画在同一张图里,横轴是 $\sigma_\alpha$ 或 $\sigma_\tau$,纵轴是定位误差。如果 RCTLS 曲线离 CRLB 在 2~3 倍以内,说明算法已经榨干了测量信息;如果差一个数量级,先去检查加权矩阵是不是没设对,而不是继续调 $\lambda$。
6.2 画误差椭圆比只看 RMSE 更能暴露系统偏差
RMSE 是一个标量,会把不同方向的误差混在一起。定位误差往往不是各向同性的——TDOA 约束强的方向误差小,DOA 约束强的方向误差大。把 Monte Carlo 的定位结果散点画出来,对协方差矩阵做特征值分解,两个特征值对应误差椭圆的长短轴,特征向量给出椭圆朝向。如果椭圆中心不在真实目标点上,说明存在系统性偏差,这时回头检查 $d_i$ 公式或者 $\lambda$ 的偏差影响。
这也是我复现这类定位算法的习惯:先看几何再看统计,先跑单次再跑批量。论文里那些漂亮的 RMSE 曲线背后,通常都藏着几何条件、正则化参数和加权矩阵的反复折腾。RCTLS 不是银弹,但它确实把 LS 在系数噪声下无能为力的那部分误差吃掉了。希望帮到你。
本文还有配套的精品资源,点击获取