MATLAB雷达系统仿真:FMCW与脉冲多普勒全链路解析
简介:一套基于MATLAB的雷达系统仿真程序,面向雷达信号处理学习者与电子信息类专业学生,覆盖信号产生、传播散射、接收处理、目标检测、干扰对抗及性能评估全链路,可用于理解脉冲、连续波、线性调频连续波等波形设计与雷达方程指标计算。压缩包共233个文件、大小约1.33MB,其中84个m脚本承担核心仿真逻辑,83个html文档提供说明与结果展示,另有c源码、dll、mexglx、mexsol等编译文件用于扩展加速,README与txt文件辅助上手。已有188人学习下载。包内实现匹配滤波、脉冲压缩、动目标检测、频率与空间分集等多种算法,并可通过Simulink构建系统级模型,便于从模块仿真延伸至硬件在环测试,为快速原型验证提供参考路径。目录结构清晰,适合直接运行、修改参数并观察结果,可作为课程设计、毕设研究或雷达方向入门进阶的操作平台。
1. 一份 MATLAB 雷达系统仿真程序:从信号生成到目标检测的完整链路
手上这份 MATLAB 雷达系统仿真程序,跑通之后可以直接输出目标的距离、速度和距离-多普勒图,链路是闭环的:波形生成、回波模拟、混频、距离维 FFT、慢时间多普勒 FFT、CFAR 检测,一个都不少。拆这类程序,最值钱的地方不是某个算法写得有多花哨,而是它把 FMCW 与脉冲多普勒两种体制里那些倍数关系(时延、差频、多普勒频率)全部落成了可运行的代码。适合两类人:刚接触雷达信号处理、想在 MATLAB 里把测距测速流程完整跑一遍的初学者,以及需要快速改波形参数做算法验证的工程人员。下面按我实际拆包的顺序来写。
2. 雷达系统仿真的模块结构与数学模型:从波形生成到检测输出,链路里到底有哪几步
拿到一个“完整”的仿真程序,第一件事不是急着运行,而是先弄清楚每个 .m 文件在链路里的位置。雷达仿真的整体结构和软件工程里的分层很像:底层是信号模型,中间是处理算法,上层是检测判决。程序里各个脚本的命名可能不一样,但功能模块逃不出这几个。
2.1 功能模块划分:先认清单个脚本在链路里的位置
我在拆这份程序时,习惯先在纸上画一条处理链,再把源码文件按功能塞进去。正常的 FMCW 雷达仿真至少有六个模块,我整理成下表方便对照:
| 功能模块 | 主要作用 | 典型输出 |
|---|---|---|
| 波形生成 | 产生 FMCW 扫频或脉冲串的复包络 | 发射信号 tx |
| 回波模拟 | 按目标距离和速度构造时延、频移 | 接收信号 rx |
| 混频处理 | 发射与接收共轭相乘得到差频信号 | 中频 IF 信号 |
| 距离维 FFT | 对单个 chirp 做快时间维傅里叶变换 | 距离谱 |
| 慢时间 FFT | 对多个 chirp 的同一距离门做傅里叶变换 | RD 谱 |
| 恒虚警检测 | 用自适应门限从频谱中提取目标峰 | 目标距离/速度列表 |
这个划分对应到代码调试非常有用。比如距离测不准,问题大概率出在回波模拟和混频两个模块;速度测不准,问题多半在慢时间 FFT 的参数上。程序中每个模块通常只依赖前一个模块的输出,所以可以分段运行、分段打印中间结果,不必每次都从头跑一遍。
有一种常见误用是拿到程序直接运行主脚本,看到图出来就认为完事了。实际上主脚本只是一条主线,各功能模块被拆成了独立函数。调试时要学会“断点定位”:先跑波形生成,确认时域波形是扫频;再跑混频,确认差频是常数;最后做 FFT,确认峰的位置。每一步都对得上,链路才成立。
2.2 回波模型:雷达方程、时延与多普勒在代码里怎么落地
雷达仿真里最核心的数学关系其实只有三个。第一个是时延,目标距离 R 对应的电磁波往返时延 τ = 2R/c,这个 2 是双程路径带来的,系数丢了整个仿真就偏了。第二个是多普勒频移,径向速度 v 对应的频移 f_d = 2v/λ,这里的 2 是电磁波反射产生的二次频移,和时延那个系数不是一回事,但经常有人搞混。
第三个是差频公式:FMCW 发射信号经过时延 τ 后与当前发射信号混频,得到的中频频率近似为 f_b = S·τ,其中 S = B/T_c 是调频斜率。把 τ = 2R/c 代进去,得到 f_b = 2BR/(c·T_c),这就是从差频频率反推距离的基本公式。这份程序里几乎所有的距离换算都在这条公式上,改参数时必须保持 B、T_c、c 的单位一致,否则误差会直接进入距离读数。
实际仿真还有一个工程约定值得注意:发射信号不必真的生成 24GHz 的载波波形,而是用复指数包络去等效。程序里通常直接生成中频或基带信号,把载频作为系数体现。这么做不损失正确性,因为混频时载波项相减约掉了,只有差频和相位项被保留。如果你看到代码里 fc 只是被当作 lambda 的分母使用,别奇怪,这是教学仿真和工程仿真的常见做法。
2.3 为什么选 MATLAB 而不是 Python:矩阵化、FFT 与工具箱的取舍
用 MATLAB 做雷达仿真,最重要的理由是复数矩阵和 FFT 的原生支持。雷达信号处理天然是二维复矩阵操作,距离维 FFT 是对矩阵的行做变换,多普勒维 FFT 是对列做变换,MATLAB 里就是一条 fft 命令加一个维度参数的事,不需要像 Python 那样先 import numpy,再考虑数据形状对齐。
另一个因素是交互式调试体验。雷达仿真最常做的事是反复改波形参数、看频谱变化,MATLAB 的变量工作区可以直接查看每个中间矩阵的尺寸和数值,配合绘图函数能快速定位异常。Python 也能做到,但需要自己搭环境,对于以算法验证为目的的资源,MATLAB 的开箱即用优势很明显。
当然,如果所在环境没有信号处理工具箱,也能跑通这份程序的核心部分,因为 FMCW 与脉冲多普勒的基础处理只需要 fft、conj、exp 这些基础函数。这也侧面说明这份程序的依赖很薄,适合在没有额外许可的环境下运行。
3. FMCW 测距测速的 matlab 核心实现:一个扫频周期的距离 FFT 与慢时间多普勒 FFT
3.1 波形参数先定这五个数
我建议在打开代码编辑器之前,先按下面的表格把波形参数定好。参数之间是互相约束的,单独调一个数很可能引发连锁问题。
| 参数 | 符号 | 示例值 | 约束来源 |
|---|---|---|---|
| 载频 | fc | 24 GHz | 决定波长,进而决定多普勒灵敏度 |
| 扫频带宽 | B | 150 MHz | 决定距离分辨率,ΔR = c/(2B) |
| 扫频周期 | Tc | 50 us | 决定最大多普勒频率,进而决定最大不模糊速度 |
| 采样率 | fs | 5 MHz | 必须大于两倍最大差频 |
| chirp 数量 | Nchirp | 128 | 决定速度分辨率和积累增益 |
一个经常被误解的点是采样率并不是越高越好。距离分辨率只由带宽 B 决定,把 fs 从 5MHz 提到 20MHz,距离分辨率不会变好,只会让每个 chirp 的采样点数变多,FFT 计算量成倍增长。采样率的下限则来自最大差频:如果最大探测距离是 R_max,最大差频 f_max = 2BR_max/(c·T_c),采样率必须大于 2f_max。这个关系在代码里体现为对 fs 的约束,新手最容易只关注带宽而忽略采样率,导致折叠。
3.2 发射、回波与混频的三段核心代码
下面这段代码是 FMCW 单目标仿真的骨架,按发射、回波、混频三部分组织,可以直接在 MATLAB 里运行:
% 波形参数 fc = 24e9; % 载频 24 GHz B = 150e6; % 扫频带宽 150 MHz Tc = 50e-6; % 扫频周期 50 us fs = 5e6; % 采样率 5 MHz c = 3e8; % 光速 lambda = c / fc; % 波长 12.5 mm S = B / Tc; % 调频斜率 3e12 Hz/s % 目标参数 R0 = 60; % 目标距离 60 m v0 = 12; % 目标径向速度 12 m/s % 快时间采样序列 Ns = round(fs * Tc); % 每个 chirp 采样点数 = 250 t = (0 : Ns-1) / fs; % 快时间轴 % 发射信号:复指数形式的线性调频 tx = exp(1j * 2 * pi * (fc*t + 0.5 * S * t.^2)); % 回波信号:时延 tau0,多普勒频移叠加在相位上 tau0 = 2 * R0 / c; % 固定时延 fd = 2 * v0 / lambda; % 多普勒频移 rx = exp(1j * 2 * pi * (fc*(t - tau0) + 0.5 * S * (t - tau0).^2 ... + fd*t)); % 回波 % 混频:发射信号共轭乘以回波,得到中频信号 if_sig = tx .* conj(rx);这里有个关键约定:程序不直接以 24GHz 载频采样,而是把载频放在复指数里参与运算,混频后自然消掉。rx 表达式中把多普勒频移写成 fd·t,是在慢时间维上做速度检测的伏笔。实际雷达里 fd 会在每个 chirp 之间产生相位变化,这里用近似形式先把单 chirp 链路跑通。
回波生成时延 tau0 用的是固定值。如果目标速度很快,理论上时延会随慢时间变化,但单个 chirp 内目标位移远小于距离分辨率,所以先做固定时延、再在慢时间维叠加多普勒,这是 FMCW 仿真里最常见且够用的近似。把 fd 写进 rx 而不是写在混频之后,是因为混频后的中频信号里频率项会包含差频与多普勒的叠加,分开处理更直观。
3.3 距离 FFT 与慢时间多普勒 FFT 的换算
混频得到的中频信号频率包含两个成分:距离差频 S·tau0 和多普勒频移 fd。要做距离和速度解耦,需要两步 FFT。第一步是距离 FFT,对单个 chirp 的中频信号做傅里叶变换,得到距离谱:
% 距离维 FFT Nfft = 2^nextpow2(Ns); % 补零到 256 点 spec = fft(if_sig, Nfft); f_axis = (0 : Nfft-1) * fs / Nfft; % 频率轴,单位 Hz % 取正频部分,找峰值 [~, k] = max(abs(spec(1 : Nfft/2))); f_b = f_axis(k); % 峰值处对应的差频频率 % 由差频反推距离 R_est = c * f_b / (2 * S); fprintf('估计距离:%.2f m\n', R_est);这段代码中,f_b 是差频频率,R_est 用 f_b = 2BR/(cT_c) 反推。注意这里还没有把 fd 分离出去,如果目标速度较快,f_b 里混有 fd 成分,距离读数会有一个偏移。在锯齿波 FMCW 里,距离和速度不能靠单 chirp 完全解耦,需要多 chirp 处理。
第二步是慢时间 FFT。把 Nchirp 个 chirp 的距离谱峰值相位按顺序排列,对慢时间维做 FFT,峰值对应的频率就是多普勒频移 fd,再换算成速度:
% 慢时间维 FFT:多普勒维 Nchirp = 128; fd_axis = (-Nchirp/2 : Nchirp/2-1) * (1/Tc) / Nchirp; % 多普勒频率轴 v_axis = fd_axis * lambda / 2; % 速度轴 % 对慢时间维做 FFT 并中心化 doppler_spectrum = fftshift(fft(range_peaks, Nchirp)); % 找峰值对应的多普勒频率 [~, kd] = max(abs(doppler_spectrum)); fd_est = fd_axis(kd); v_est = fd_est * lambda / 2; fprintf('估计速度:%.2f m/s\n', v_est);range_peaks 是之前 Nchirp 个 chirp 距离谱峰值处的复数值,这一步就是完整的慢时间处理。速度换算公式 v = fd·λ/2 来自多普勒频移定义,注意这里的 2 是反射系数,不是双程时延那个。
4. 脉冲多普勒模式与距离-多普勒图:把快时间矩阵转成 RD 谱的两步 FFT
FMCW 只是这份程序支持的一种波形。程序里还带脉冲多普勒模式,对应的核心产物是距离-多普勒图(RD 谱)。RD 谱是雷达信号处理里最具信息量的一张图:横轴是距离,纵轴是速度,每个亮点代表一个目标。理解这张图的关键是把数据维度弄清楚。
4.1 快时间与慢时间:一个 chirp 和一帧的关系
雷达接收数据在时间上有两个维度。快时间是单个 chirp 内部的采样序列,持续时间只有 50us,采样率 5MHz 下是 250 个点;慢时间是不同 chirp 之间的索引,128 个 chirp 构成一帧,帧时间 T_frame = Nchirp·Tc = 6.4ms。快时间维蕴含距离信息,慢时间维蕴含多普勒信息。
这两个维度的角色容易混淆。做距离 FFT 时是对快时间维做变换,得到每个 chirp 的距离谱;做多普勒 FFT 时是对慢时间维做变换,同一距离单元在不同 chirp 之间的相位变化被提取出来。程序里这一步通常组织成矩阵:行是慢时间索引,列是快时间 FFT 结果,行列顺序不同直接影响后面 RD 谱的方向。
4.2 从快时间距离谱到 RD 谱的两步 FFT
组织 RD 谱的代码模式很固定。先用循环把每个 chirp 的距离 FFT 结果存入矩阵,再对矩阵的慢时间维做 FFT,最后做 fftshift 把零频移到中心:
% 生成 Nchirp x Nfft/2 的距离谱矩阵 range_profile = zeros(Nchirp, Nfft/2); for m = 1 : Nchirp % 每个 chirp 的发射和回波,目标速度在慢时间维体现为相位旋转 tx_m = exp(1j * 2 * pi * (fc*t + 0.5 * S * t.^2)); tau_m = 2 * (R0 + v0 * m * Tc) / c; % 慢时间时延变化 rx_m = exp(1j * 2 * pi * (fc*(t - tau_m) + 0.5 * S * (t - tau_m).^2)); if_m = tx_m .* conj(rx_m); % 混频 range_profile(m, :) = fft(if_m, Nfft); % 取正频部分存入矩阵 end % 慢时间维 FFT,得到距离-多普勒图 rd_map = fftshift(fft(range_profile, Nchirp, 1), 1); % 频率轴和距离轴换算 range_axis = (0 : Nfft/2-1) * fs / Nfft * c / (2 * S); fd_axis = (-Nchirp/2 : Nchirp/2-1) / Nchirp * (1 / Tc); v_axis = fd_axis * lambda / 2;这段代码里最需要注意的有两点。第一,慢时间维变化体现在 tau_m 随 m 增大,目标靠近时延减小,回波相位相对每个 chirp 产生线性旋转。第二,fft 的第三个参数 1 指定沿第一维(慢时间)做变换,很多人在这一步把维度写错,导致 RD 谱看起来像噪声。
tau_m 里加了 v0·m·Tc,这相当于把目标的径向运动分解到每个 chirp 的时延上。对慢时间维做 FFT 后,这个线性变化的相位会形成一个谱峰,位置对应多普勒频移。把时延变化和相位变化放在一起,是 FMCW 仿真里比较贴近物理实际的建模方式,比单纯在单 chirp 里叠加 fd 更利于理解 RD 谱的成因。
4.3 从 RD 谱读目标:峰值坐标、分辨率和边界
RD 谱生成后,读目标的方法是先找二维峰值,再把峰值坐标映射到距离轴和速度轴。峰值所在行对应多普勒频率,也就是速度;所在列对应差频,也就是距离。这两个换算在上面的 range_axis 和 v_axis 里已经完成,直接索引即可。
分辨率的边界由波形参数决定。距离分辨率 ΔR = c/(2B),带宽 150MHz 时约为 1m,两个目标距离差小于 1m 时在距离维上会合成一个峰。速度分辨率 Δv = λ/(2·T_frame),帧时间 6.4ms 时约为 0.98m/s,约等于 1m/s。注意帧时间越长速度分辨率越好,但帧时间变长会让数据量变大、目标在帧内的加速度影响积累效果,这是系统设计上的取舍。
RD 谱还有一个不可回避的边界:最大不模糊速度。多普勒频率超过 PRF/2 时会发生速度回卷,表现为速度符号反转或折叠到另一侧。所以改 PRF 时,既要看最大探测距离是否满足,也要看最大速度是否落入无模糊区间。这份程序里这两个约束都写在参数区注释里,属于验证参数合理性时最先要检查的项目。
5. 波形参数排查:采样率、FFT 点数与速度模糊的三个高频翻车点
5.1 距离读数整体偏移:双程时延的系数丢了
现象:仿真里设置目标距离 60m,输出 29.85m 或 119.4m,读数与真值差一个固定比例。解决办法是检查代码里时延和距离换算公式中的 2。
原因:时延 τ = 2R/c 里的 2 是电磁波往返路程,距离换算 R = c·f_b/(2S) 里的 2 对应差频公式推导中的系数。这两个地方漏掉或写错,读数不是偏小一半就是偏大一倍。
解决:先打印混频后的差频频率,用手算的 f_b = 2BR/(c·T_c) 对比。如果仿真频率对而距离不对,问题出在换算公式;如果频率本身不对,问题出在时延构建。我给每个涉及 2 的公式单独写注释,避免下次再犯。
5.2 目标峰直接消失:采样率低于最大差频
现象:目标距离超过某个值后,距离谱上找不到峰,或者峰出现在错误的位置,像是被“截断”了。
原因:最大差频 f_max = 2BR_max/(c·T_c) 超过采样率的一半,中频信号发生频谱混叠。混叠后峰位不再是真实的差频,距离解算必然出错。
解决:把 fs 提高到至少两倍最大差频。经验做法是先算理论最大差频,再留 1.5 倍余量。注意提高 fs 会增加每 chirp 点数,FFT 计算量上升,但这是保证正确性的必要代价。
5.3 速度正负读反或回卷:多普勒频率超过 PRF/2
现象:目标径向速度 12m/s,RD 谱显示的却是接近 -12m/s 的值,或者速度读数明显不对。
原因:慢时间采样率就是 chirp 重复频率 PRF = 1/Tc。多普勒频移超过 PRF/2 时发生欠采样,速度折叠到负方向或无模糊区间之外。
解决:按最大目标速度反推 Tc 上限:Tc ≤ λ/(4·v_max)。如果速度范围很大,需要减小 Tc,这会压缩最大探测距离,所以工程上常需要牺牲距离或使用多重 PRF。程序里如果保留了 PRF 参数,改完 Tc 后要重新核对最大距离是否仍满足要求。
5.4 FFT 补零与复数混频:两个“看起来对但结果错”的细节
现象:对中频信号补了很多零,频谱峰看起来变细了,但目标距离精度并没有提升;或者复数混频后频谱左右出现对称的两个峰,不知道该取哪个。
原因:补零只是对频谱做插值,让峰值位置更精细,并不能提高真实分辨率。两个对称峰多半是实信号混频后取实部丢掉了负频率信息,或者 I/Q 通道幅度不一致导致镜像残留。
解决:用复信号混频并保留正频部分,峰值位置用抛物线插值或补零提升精度。镜像问题要检查 I/Q 通道的幅度和相位是否平衡,程序里如果使用单通道实信号建模,就要做好只取单边谱的准备。
6. 进阶技巧:CA-CFAR 自动门限与注入目标的仿真自检
6.1 CA-CFAR:把固定门限换成自适应门限
距离和速度都能测了,下一步是让程序自己找目标。仿真里最常见的做法是设置一个固定幅度门限,目标明显时没问题,但换一组参数后门限可能失效。CA-CFAR(单元平均恒虚警)用待检测单元两侧的参考单元平均值作为底噪估计,再乘以系数得到门限,能自适应当前噪声水平,是雷达检测的标准做法:
function th = ca_cfar(x, ng, nr, pfa) % x 输入频谱,ng 单侧保护单元数,nr 单侧参考单元数,pfa 虚警概率 Nref = 2 * nr; % 参考单元总数 alpha = Nref * (pfa^(-1/Nref) - 1); % CFAR 系数 th = zeros(size(x)); for i = ng + nr + 1 : length(x) - ng - nr ref = [x(i-ng-nr : i-ng-1), x(i+ng+1 : i+ng+nr)]; % 两侧参考窗 th(i) = alpha * mean(abs(ref)); % 门限 = 系数 * 平均噪声功率 end end保护单元是必要的:目标峰展宽后可能漏进参考窗,把门限抬得过高,导致真实目标被过滤。ng 通常取 2 到 4 个单元,nr 取 8 到 16 个单元;pfa 取 1e-4 或 1e-6,越小门限越高。这个函数可以直接替换固定门限,放在距离 FFT 和慢时间 FFT 之后,对 RD 谱的每一行做处理。
6.2 注入已知目标做回归验证:养成仿真自检的习惯
CFAR 跑通之后,我最推荐做的是注入已知目标验证。做法是设置一组 R_test 和 v_test,用程序生成多目标回波,经过完整链路后把估计值与真值对比。这相当于给仿真做回归测试,改过参数后跑一遍就知道有没有引入偏差:
R_test = [30 50 70]; v_test = [5 -8 12]; % 按这组真值生成多目标回波,跑完整链路 % 检查每个目标输出:距离误差应小于一个距离分辨率,速度误差应小于一个速度分辨率拿着这张误差表去检查参数是最直接的。误差大于分辨率说明链路里某个近似失效,常见原因包括目标间距小于一个分辨单元、多普勒超出无模糊区间、或采样率贴近混叠边界。从那以后我每次改完波形参数,都强制把“最大距离、最大速度、采样率、FFT 点数”四个数写在一张纸上对一遍,再跑仿真确认误差表,不会再凭感觉调参数。希望这份拆解能帮你在自己的项目里少走几步弯路。
本文还有配套的精品资源,点击获取