多普勒模糊原理与Matlab仿真：脉冲雷达测速中的频谱折叠与解模糊方法-平芜编程栈

做雷达信号处理的人，十有八九都遇到过这种场景：明明仿真里的目标速度已经到几十米每秒了，多普勒谱上却在一个很低的频率位置冒出一个峰值，看起来像是一个“慢速目标”。我最早做脉冲多普勒雷达实验时也在这个问题上栽过跟头，排查半天发现不是代码写错了，而是多普勒模糊在捣乱。

多普勒模糊是雷达、声呐、超声成像里绕不开的经典问题。这个仿真实验的核心，就是用Matlab把多普勒模糊现象完整复现出来，通过构造脉冲串回波、做慢时间维FFT，直观看到真实多普勒频率如何被“折叠”到低频区间，再通过解模糊处理还原目标真实速度。内容覆盖脉冲多普勒雷达的信号模型、欠采样原理、模糊速度推导、谱峰检测和速度解算全套流程。适合正在学雷达信号处理的本科生、研究生，也适合做通信或声学仿真的人参考，代码可以直接在Matlab里跑起来，改改参数就能当模板用。

1. 多普勒模糊是什么：先把这个现象讲透

1.1 雷达测速为什么躲不开多普勒效应

雷达测速的原理其实一句话就能说清：发射电磁波照到目标上，目标有径向运动时，回波频率会发生偏移，这个偏移量就叫多普勒频率。对一个径向速度为v的目标，多普勒频率为：

fd = 2v / λ

其中λ是雷达工作波长。从这个式子能看出两件事：一是速度越快，多普勒频率越高；二是波长越短，相同速度对应的多普勒频率越高。比如我后文仿真里用10GHz载频，波长0.03m，一个19.5m/s的目标，多普勒频率就是1300Hz。这个频率放到连续波雷达里很好测，直接做频谱分析就能看到。

但脉冲雷达不是这样工作的。脉冲雷达发射的是一个个离散脉冲，接收目标回波也是在脉冲间歇里采样，本质上是在用一段脉冲串去观测目标的相位变化。那么问题来了：这种脉冲化的观测方式，对多普勒频率的测量范围有没有限制？

有，而且限制非常严格。

1.2 慢时间采样与频谱折叠

雷达里有两个时间维度的概念，新手很容易搞混。快时间指一个脉冲内部的采样时间，用来分辨距离；慢时间指脉冲序号对应的时间，也就是每个脉冲到达目标的时刻序列。目标运动导致的多普勒信息恰恰藏在慢时间维的相位变化里。

在慢时间维上，雷达的采样率就是脉冲重复频率PRF。按奈奎斯特采样定理，信号不模糊可观测的频带范围是[-PRF/2, PRF/2]。目标真实多普勒频率一旦超过这个范围，就会发生频谱折叠，像是把一个更高的频率搬移到了低频区间。这种现象在周期信号采样里非常典型，和音频信号采样率不足导致混叠是一个道理，只不过在雷达里影响的是速度测量。

具体折叠公式可以写成：

fd_obs = mod(fd_true + PRF/2, PRF) - PRF/2

也就是说，无论真实多普勒有多高，我们最终看到的表观多普勒频率永远落在[-PRF/2, PRF/2]之间。后文仿真中真实多普勒是1300Hz，PRF是1000Hz，代入公式得到表观频率300Hz，正好是1300Hz减去一个1000Hz的结果。

1.3 最大不模糊速度的推导

把折叠公式和波长关系联立起来，就能推出一个重要指标：最大不模糊速度。当目标的真实多普勒频率等于PRF/2时，恰好是观测范围的边界，超过这个边界就开始折叠，所以最大不模糊速度为：

v_max = λ·PRF / 4

以我仿真的参数计算，λ=0.03m，PRF=1000Hz，v_max只有7.5m/s。而目标速度是19.5m/s，超过了v_max约2.6倍，所以必然产生模糊。

这个公式还揭示了一个雷达设计中非常经典的矛盾：想测高速目标，就要提高PRF，但提高PRF后脉冲间隔变短，最大不模糊距离c/(2PRF)又会缩小。距离测量和速度测量在PRF这个参数上此消彼长，怎么取舍是雷达总体设计的大课题，而这个仿真实验恰好能把矛盾中的一半——速度模糊——完整呈现出来。

2. 仿真实验设计：参数怎么定、信号怎么建模

2.1 场景设定与关键参数选择

仿真实验的第一步是定参数。这一步看着简单，实则很考验对原理的理解，参数选得不好，要么模糊现象不明显，要么代码跑起来数据量太大。

我选的载频是10GHz，对应波长0.03m。目标初始斜距3000m，径向速度19.5m/s，理论多普勒频率1300Hz。PRF设为1000Hz，这样观测频带是[-500Hz, 500Hz]，能保证目标多普勒处于模糊状态。脉冲积累数取64个，这个数量在速度和频率分辨率之间比较均衡，再多仿真矩阵会变大，再少频谱峰值不够锐利。

快时间采样率设为10MHz，采样512个点，对应的距离门宽度为c/(2Fs)=15m，距离观测窗口覆盖7680m，目标放中间没问题。这里有个细节值得注意：在64个脉冲积累时间内，目标实际只移动了约1.25m，不到一个距离门，所以目标始终落在同一个距离门上。这样处理起来非常方便，后文直接取目标距离门的慢时间信号做FFT即可，完全不用担心距离徙动破坏多普勒谱。

2.2 回波建模：从公式到矩阵

回波建模是本实验的核心环节。点目标回波可以用一个简明的表达式描述：快时间维上是发射脉冲的延迟复制，慢时间维上附加一个多普勒相位旋转。在基带等效模型里，回波延时为：

τ(m) = 2(R0 + v·m·PRT) / c

其中m是慢时间脉冲序号。目标回波的复振幅在快时间维上位于τ(m)对应的距离门，相位在慢时间维上按exp(j2π·fd·m·PRT)旋转。

把这个公式转成代码思路就非常直接了：构造一个Npulse×Nfast的复数矩阵，每一行是一个脉冲的快时间采样，每一列是一个距离门在慢时间维上的时间序列。目标回波就是在目标所在距离门的位置放入一个复数值，该复数的幅角随时间按多普勒频率旋转。

这种建模方式相比逐采样点生成发射波形再求卷积，要高效得多，而且物理含义清晰：距离维由快时间采样决定，速度维由慢时间相位变化决定，两者天然正交，便于后续做二维FFT处理。

2.3 为什么用快时间-慢时间数据矩阵

刚接触雷达仿真的朋友可能不太理解，为什么要费劲构造这样一个二维矩阵，直接把回波当成一维信号处理不行吗？

不行的根本原因在于脉冲雷达天生就是二维采样结构。一个脉冲对应一个时间点，这个时间点内目标回波包含了距离信息；多个脉冲串起来，同一距离门在不同脉冲间的相位变化才包含速度信息。如果把回波简单拉成一维序列，快时间和慢时间会混在一起，后续的距离-多普勒二维处理就无从谈起。

用快时间-慢时间矩阵还有一个好处：可视化直观。把矩阵做二维FFT之后，横轴是距离、纵轴是多普勒频率，目标和杂波在二维平面上呈现为明亮的峰值，距离和速度同时读出。这种处理方式就是教科书里常说的MTD（动目标检测）的基本结构，掌握这个矩阵思维，后面学空时自适应处理都会顺畅很多。

3. Matlab仿真实现：从零搭一个可运行实验

3.1 参数初始化代码

下面这段代码是整个实验的骨架。建议直接在Matlab里新建脚本，按顺序粘贴运行，R2016之后的版本基本都能跑通。

clear; close all; clc; %% 系统参数 fc = 10e9; % 载频 10 GHz c = 3e8; % 光速 lambda = c / fc; % 波长 0.03 m PRF = 1000; % 脉冲重复频率 1000 Hz PRT = 1 / PRF; % 脉冲重复周期 1 ms Npulse = 64; % 慢时间积累脉冲数 Fs = 10e6; % 快时间采样率 10 MHz Nfast = 512; % 快时间采样点数（距离门数） %% 目标参数 R0 = 3000; % 初始斜距 3 km v = 19.5; % 径向速度 m/s fd_true = 2 * v / lambda; % 真实多普勒频率 1300 Hz %% 快时间轴与距离门 t_fast = (0:Nfast-1) / Fs; range_axis = c * t_fast / 2; %% 慢时间轴 t_slow = (0:Npulse-1) * PRT;

参数注释里的单位我全标清楚了，这是好习惯，尤其是涉及频率和速度换算时，单位不统一很容易出低级错误。

3.2 构造回波与加噪

回波构造按前文说的思路实现。目标时延随时间变化体现在距离门索引的微小移动上，而每个脉冲上目标的复振幅按多普勒频率旋转：

%% 每个慢时间脉冲对应的目标距离门 tau = 2 * (R0 + v * t_slow) / c; delay_idx = round(tau * Fs) + 1; %% 构造快时间-慢时间回波矩阵 s = zeros(Npulse, Nfast); for m = 1:Npulse s(m, delay_idx(m)) = exp(1j * 2 * pi * fd_true * t_slow(m)); end %% 加复高斯白噪声，目标距离门处信噪比 20 dB SNR_dB = 20; noise_power = 10^(-SNR_dB/10); noise = sqrt(noise_power/2) * (randn(Npulse, Nfast) + 1j*randn(Npulse, Nfast)); s = s + noise;

这里我刻意把延迟索引的移动幅度控制在一个距离门以内。验算一下：64个脉冲积累时间0.064s，目标移动1.25m，对应时延变化约8.3ns，换算成距离门只有0.083个，于是delay_idx在64个脉冲里几乎不变。这在代码里是有意为之，核心目的是让目标能量集中在一个距离门内，多普勒频率的模糊特性就能干净地展示出来。如果目标速度大到积累时间内跨多个距离门，就得额外做距离徙动校正，那就超出这个实验的范畴了。

3.3 多普勒FFT与频谱绘图

对慢时间维做FFT是观察多普勒模糊的关键一步。注意FFT的方向是沿着矩阵的第一维，也就是脉冲序号方向，这对应的是慢时间。做完之后用fftshift把零频挪到频谱中心，频率轴也同步平移：

%% 沿慢时间维做FFT，得到距离-多普勒谱 S = fftshift(fft(s, Npulse, 1), 1); %% 多普勒频率轴 fd_axis = (-Npulse/2 : Npulse/2-1) / Npulse * PRF; %% 取出目标所在距离门的多普勒谱 [~, idx_r] = max(abs(s(1, :))); figure; plot(fd_axis, 20*log10(abs(S(:, idx_r)) / max(abs(S(:, idx_r)))), 'LineWidth', 1.2); xlabel('多普勒频率 (Hz)'); ylabel('归一化幅度 (dB)'); title('目标距离门的多普勒频谱'); grid on;

运行后频谱峰值会出现在300Hz附近，而不是真实值1300Hz。这个300Hz就是慢时间欠采样导致的折叠频率。我第一次跑出这个结果时还有点怀疑代码写错了，后来把观测范围画出来才意识到，1300Hz落在1000Hz采样率下必然折到300Hz，物理过程完全正确。

再画一张距离-多普勒二维图，效果更直观：

figure; imagesc(range_axis, fd_axis, 20*log10(abs(S))); xlabel('距离 (m)'); ylabel('多普勒频率 (Hz)'); title('距离-多普勒二维谱'); colorbar;

二维图上能清楚看到目标峰值在(3000m, 300Hz)附近。距离坐标基本和设定的R0一致，多普勒坐标则把模糊现象摆在了面前。

4. 解模糊处理：从模糊频率还原真实速度

4.1 目标检测与峰值提取

看到模糊频率之后，接下来的问题自然是怎么“解模糊”，也就是从表观多普勒频率反推出真实速度。第一步先是峰值检测。

严格工程化的做法是用CFAR检测器，在距离-多普勒二维平面上自适应求阈值，把超过阈值的峰值位置提取出来。教学实验不用那么复杂，直接在目标所在距离门的频谱里找最大峰值就行：

%% 在目标距离门的多普勒谱中找峰值 [S_max, idx_max] = max(abs(S(:, idx_r))); fd_obs = fd_axis(idx_max); fprintf('检测到的表观多普勒频率: %.2f Hz\n', fd_obs);

如果是低信噪比场景，可以改用全二维平面找全局最大值，或者先设定一个幅度阈值再找极值点。这个实验SNR有20dB，目标峰值比噪声高几个数量级，简单峰值检测完全够用。

4.2 模糊数估计与速度还原

峰值检测拿到fd_obs=300Hz之后，需要判断它对应哪个模糊周期。真实多普勒频率满足以下关系：

fd_true = fd_obs + k·PRF

其中k是整数，称为模糊数。问题是k可以取任何整数值，300Hz、1300Hz、2300Hz全都不区别，怎么确定是哪一个？

这就是解模糊的核心难点。实际雷达系统里通常有几个办法：一是利用目标跟踪的先验信息，知道目标大致速度范围，直接把不合理的k排除；二是用两个不同PRF交替测量，联合求解；三是用目标跨距离门的速率粗测速度，再判断k。前两种是系统工程上最常用的，第二种我放到下一小节单独讲。

这里先说一个教学上最简单的方法——利用距离门迁移量。回波数据里其实已经包含了目标慢速移动的信息：虽然目标在64个脉冲内没跨过完整距离门，但相位上是连续变化的。如果积累脉冲数足够多，目标跨过的距离门数就可以反算出粗速度：

%% 用累积时间内的距离门跨度粗测速度 delta_gate = delay_idx(end) - delay_idx(1); v_coarse = delta_gate * (c / (2*Fs)) / (Npulse * PRT);

把这个粗测速度和观测频率结合起来，模糊数就可以唯一确定了。比如粗测速度约20m/s，对应真实多普勒约1333Hz，那么k=1显然比k=0、k=2更合理，最终得到真实速度19.5m/s。

4.3 多PRF联合解模糊扩展实验

多PRF解模糊是工程上最实用的方案，原理深刻但并不难懂。核心思想是：用两个不同的PRF分别测量同一个目标，得到两个不同的折叠频率，再根据中国余数定理的思想恢复真实频率。

我用另一组参数验证过这个思路。设PRF1=1000Hz，观测到fd1=300Hz，那么：

fd_true = 300 + k1·1000

设PRF2=1400Hz，观测范围是[-700Hz, 700Hz]，1300Hz在这个范围之外，折叠后得到fd2=-100Hz，那么：

fd_true = -100 + k2·1400

枚举k1和k2，找同时满足两个条件的频率解。k1=1时得到1300Hz，同时k2=1也得到1300Hz，两边对上，1300Hz就是真实多普勒频率。

下面这段小代码可以自动完成枚举：

PRF_set = [1000, 1400]; fd_obs_set = [300, -100]; fd_candidates = []; for k = -5:5 fd_candidates = [fd_candidates, fd_obs_set(1) + k * PRF_set(1)]; end for k = -5:5 fd = fd_obs_set(2) + k * PRF_set(2); if any(abs(fd_candidates - fd) < 1e-6) fprintf('解模糊后的真实多普勒频率: %.2f Hz\n', fd); fprintf('目标真实速度: %.2f m/s\n', fd * lambda / 2); end end

这段代码写得比较直白，实际工程里可以用扩展的欧几里得算法直接求同余方程，不用枚举。教学实验用枚举反而更清楚，因为能直观看到两个候选集合的交点只有一个。

这个扩展实验强烈建议读者跑一下，把PRF_set和fd_obs_set换成自己的参数，能加深对模糊折叠公式的理解。

5. 常见问题与调试经验

5.1 为什么峰值位置总差一点：栅栏效应

有读者按上面的代码跑完，发现峰值多普勒频率不是正好300Hz，而是295Hz或305Hz之类，于是怀疑代码有bug。其实这大概率是FFT栅栏效应造成的。

FFT输出的频率点离散分布，频率分辨率是PRF/Npulse，也就是1000/64=15.625Hz。真实多普勒1300Hz折叠后的300Hz不一定正好落在某个离散频率点上，峰值就会被“夹”在两个FFT bin之间，看起来整体偏移几个Hz。这并不影响理解现象，但如果要做精确测速，就要用插值或频谱细化技术。最简单的方式是把积累脉冲数增加，比如改成128或256，频率分辨率提升，峰值位置就更准。

5.2 频谱泄露干扰旁瓣怎么压

慢时间FFT的点数有限，相当于给无限长的多普勒信号加了一个矩形窗，频域上会产生旁瓣。噪声不高时旁瓣无伤大雅，但如果两个目标速度相近，旁瓣可能互相掩盖，甚至产生虚假峰值。

处理方式和时域窗函数一样，在慢时间FFT之前给数据乘一个窗函数，比如汉明窗或布莱克曼窗。代价是主瓣会展宽，需要自己权衡。教学实验中如果只是单目标场景，加不加窗影响不大，我倾向于不加窗以保留原始谱的锐利形状；但一旦扩展到多目标场景，加窗几乎是必选项。

5.3 盲速与多普勒奇偶性问题

这个坑比较隐蔽，属于“不遇到一次根本想不起来”的问题。当目标真实多普勒频率恰好等于PRF的整数倍时，折叠后的表观频率是0Hz，目标看起来完全静止，雷达会彻底丢失目标。这个速度就叫盲速。

比如本节仿真里如果目标速度是75m/s，真实多普勒5000Hz，在PRF=1000Hz下折叠为0Hz，峰值出现在零频，和目标完全不运动没有任何区别。实际雷达里盲速是硬伤，解决思路大多是采用参差重复频率，让不同PRF的盲速不重合，联合观测避免死角。

建议读者在仿真里把目标速度改成75m/s试一次，你会看到多普勒谱上0Hz处出现一个峰值，而目标真实速度却远不止静止。这个现象挺震撼的，比看十页公式记忆都牢。

5.4 代码运行效率和内存优化

这个教学实验矩阵只有64×512个复数，Matlab跑起来毫无压力。但如果把Npulse提到上千、Nfast提到上万，循环构造回波的效率短板就会显现。我这里提供一个优化思路：整个回波矩阵的构造其实可以向量化，核心是将相位历史和距离门索引用矩阵运算直接生成。不过对于理解原理而言，循环优先，代码可读性更重要，不要一上来就追求向量化而牺牲可读性。

另外真实验证时注意Matlab版本差异。fftshift的用法、复数随机数生成方式不同版本略有差异，如果代码报错，优先检查这两个函数有没有按当前版本的语法调用。