简介面向雷达信号处理学习者的 MATLAB 模糊函数仿真代码包聚焦矩形脉冲、高斯脉冲与 LFM 信号三类典型波形的模糊函数分析适合通信、雷达与电子对抗方向的学生和研究人员用于课堂实验、课程设计或论文复现。资源共 6 个 m 文件压缩包约 4KB主要包含各波形函数的实现脚本及统一封装的 ambiguity() 函数可一次性输出模糊函数图、模糊度图、距离模糊函数图与速度模糊函数图另有 LFM 信号调用示例便于直接运行对照。作者将不同脉冲信号的模糊函数计算逻辑分文件解耦封装后的函数支持按需传入参数直接复现省去从零编写与调试绘图代码的时间通过对比三类波形在时延—多普勒平面的能量分布可直观理解距离分辨率、速度分辨率与旁瓣特性之间的差异。目前已有 4217 人学习适合作为模糊函数入门仿真与结果核对的参考素材。1. 一条斜脊线值不值得信MATLAB 仿真模糊函数要先看清什么雷达波形评审会上最常被追问的一句是目标回波同时带上几微秒时延和几十千赫多普勒时匹配滤波还剩多少峰值MATLAB 仿真模糊函数就是回答这个的工具——把发射波形和它自身的时延、多普勒副本做相关铺出一张 τ-fd 二维图。矩形脉冲、高斯脉冲、LFM 信号这三类基础波形在图上完全是三副面孔矩形是十字形 sinc 栅格高斯是一团没有旁瓣的椭圆包LFM 是一条斜脊。做波形论证、写课程设计、给系统挑波形的人用同一套代码就能把三者摆在一张图里比分辨率、比旁瓣、比耦合。2. 定义式落到网格上三种波形共用的一套 ambg 计算骨架2.1 连续定义、离散化与三个必须钉死的约定模糊函数描述的是波形与自身副本的二维相关程度工程上常用下面这个形式chi(tau, fd) ∫ s(t) * conj(s(t - tau)) * exp(-j*2*pi*fd*t) dt其中s(t)是复基带发射信号tau是时延fd是多普勒频移。取模平方得到的是模糊函数模值取模本身是幅度模糊函数论文里画图多数用20*log10(|chi|/max)。这个式子的物理含义很直白先把信号延迟tau再乘上一个多普勒相位旋转然后和原信号做内积。峰值落在(0,0)因为只有零时延零多普勒时才完全对齐。落到离散域有三件事必须先定定不下来后面全是玄学。第一是多普勒符号。exp(-j*2*pi*fd*t)和exp(j*2*pi*fd*t)只差一个共轭但在图上表现为整张图上下翻转。MATLAB 的fft自带exp(-j*2*pi*k*n/N)的核所以用fft实现时应该配exp(-j...)的定义这样频率栅格和fftshift后的轴天然对齐不用额外取共轭。第二是共轭顺序。s(t)*conj(s(t-tau))和conj(s(t))*s(t-tau)互为共轭模值一样但相位相反。做模值图无所谓一旦要拿相位做相干积累就必须统一。第三是时间零点。离散序列s[n]的第 0 个采样对应t0还是t-T/2会影响你画出来的tau轴是不是以零为中心。我一般把时间轴直接生成成t ((0:N-1) - floor(N/2))/fs让零时刻落在数组中间后面所有切片都不需要再搬移。离散化本身不复杂tau m/fsm取整数也就是只在采样点整数倍上取时延。有人会问能不能插值出亚采样时延可以但插值相当于给信号额外加了一级低通会把模糊函数的高频细节抹掉做旁瓣分析时反而失真。正确做法是把采样率提上去fs 8BB 为信号带宽让1/fs足够小时延栅格自然就细了。多普勒维用 DFT 实现栅格间隔是fs/Nfft。这里Nfft只是补零长度补零不增加信息只是把 DTFT 的采样点加密让图看起来连续。真正决定多普勒分辨率的是信号时宽T量级在1/T跟Nfft无关——这一点后面避坑章节还会再提一次。计算量上直接按定义做二重循环是O(N^2 * Nfft)N2000时基本跑不动。改成「固定一个时延、对整段乘积做一次 FFT」之后复杂度降到O(N * Nfft * log(Nfft))单个波形的整张图在普通笔记本上不到一秒。2.2 一个 ambg 函数吃三种波形把上面的约定写成一个独立函数输入只有信号、采样率和两个网格参数输出模糊函数矩阵和两条坐标轴。三种波形共用它不做任何分支判断。function [AF, tau, fd] ambg(s, fs, M, Nfft) % ambg 离散模糊函数 % chi(tau,fd) sum_n s[n]*conj(s[n-m])*exp(-j*2*pi*fd*n/fs) % s : 复基带采样列向量长度 N % fs : 采样率 Hz % M : 单边时延点数tau (-M:M)/fs要求 M N-1 % Nfft : 多普勒维 FFT 点数建议 8*N s s(:); N numel(s); m -M:M; tau m(:). / fs; % 1 x (2M1) fd (-Nfft/2 : Nfft/2-1). * (fs/Nfft); % Nfft x 1升序 AF zeros(Nfft, numel(m)); for k 1:numel(m) mk m(k); ss zeros(N,1); % 有效重叠区间n 从 max(1,mk1) 到 min(N,Nmk) n1 max(1, mk1); n2 min(N, Nmk); idx (n1:n2).; ss(idx) s(idx) .* conj(s(idx - mk)); % s[n]*conj(s[n-m]) AF(:,k) fftshift(fft(ss, Nfft)); % 多普勒维一次 FFT end endss的长度固定为N只在重叠区间上填值其余补零。这一步决定了模糊函数沿tau方向自然衰减——时延越大重叠越少能量越低矩形脉冲画出来就是标准的三角形包络。如果图省事用xcorr拿到2N-1长度的结果再截断很容易把零时延对到错误的索引上这是新手最常翻的车。fftshift必须在fft之后、赋值之前做这样AF的行顺序和fd的升序严格对应后面imagesc和切片索引都不用再换算。下面生成三种波形。参数取T20 us、B5 MHz、fs8B40 MHz时宽带宽积 100属于典型的窄带脉冲压缩场景。T 20e-6; B 5e6; fs 8*B; N round(T*fs); t ((0:N-1). - floor(N/2)) / fs; % 零时刻居中 K B/T; % 调频率 Hz/s % 1) 矩形脉冲 s_rect double(abs(t) T/2); % 2) 高斯脉冲截断到 ±4 sigma sigma T/8; a 1/(2*sigma^2); s_gauss exp(-a*t.^2); % 3) LFM复基带 s_lfm exp(1j*pi*K*t.^2) .* double(abs(t) T/2); M N-1; Nfft 2^nextpow2(8*N); [AF1, tau, fd] ambg(s_rect , fs, M, Nfft); [AF2, ~ , ~ ] ambg(s_gauss, fs, M, Nfft); [AF3, ~ , ~ ] ambg(s_lfm , fs, M, Nfft);M N-1把时延铺满整个脉宽|tau| T之后没有重叠模糊函数理论上为零图的两侧会自然收边。Nfft取到8N是为了让多普勒栅格足够密fd间隔约fs/(8N) 1/(8T)比理论分辨率1/T细八倍画出来的脊线才平滑。2.3 坐标轴、归一化与出图算完之后最容易出问题的是显示环节。模糊函数的动态范围通常在 60 dB 以上线性刻度画出来只有中心一个亮点什么都看不见直接imagesc(abs(AF))还会因为imagesc默认 y 轴向下而在视觉上把图上下颠倒。function show_af(AF, tau, fd, T, name) AFn 20*log10(abs(AF)/max(abs(AF(:))) eps); imagesc(tau*1e6, fd/1e6, AFn); axis xy; % y 轴向上fd 递增 colormap(jet); caxis([-60 0]); colorbar; xlabel(\tau (\mus)); ylabel(f_d (MHz)); title(name); xlim([-1.2*T 1.2*T]*1e6); ylim([-4/T 4/T]/1e6); % 只看有意义的多普勒范围 endeps是为了防止某个点恰好为零时log10(0)变成-Inf那样整张图的caxis会被拉爆最后满屏噪点。caxis([-60 0])固定动态范围三个信号才能横向对比旁瓣高低如果每张图各自动态缩放看起来高斯和 LFM 的旁瓣差不多结论就完全错了。ylim限制在±4/T也有讲究。多普勒轴的有效范围其实只到±1/T量级再往外全是零深区铺满整个±fs/2只会让主瓣挤成一条线。截到四倍分辨率带宽主瓣和第一旁瓣的比例关系才看得清。3. 矩形脉冲与高斯脉冲先用解析式把曲线钉死再谈图3.1 矩形脉冲三角形自相关叠上 sinc矩形脉冲的模糊函数有闭式解|tau| T时|chi(tau, fd)| (T - |tau|) * |sinc( (T - |tau|) * fd )|其中sinc(x) sin(pi*x)/(pi*x)。这个式子把两件事叠在了一起前面那个(T-|tau|)是三角形的自相关包络后面那个 sinc 是对每个时延切片沿多普勒方向的调制。tau0时退化为T*|sinc(T*fd)|就是矩形脉冲频谱的模fd0时退化为T-|tau|就是标准三角形。理解这个结构对判图很关键。矩形脉冲的模糊函数图是一个沿两个轴都展开的十字栅格中心峰周围有一圈一圈的旁瓣。这些旁瓣不是数值误差是矩形硬截断带来的频谱泄漏在二维上的体现。第一旁瓣在fd ≈ ±1.5/T处电平约 -13.2 dB跟匹配滤波的距离旁瓣是同一个来源。% 零多普勒切片与解析式对比 R_sim abs(AF1(:, M1)); R_sim R_sim / max(R_sim); R_th max(0, 1 - abs(tau)/T); % 三角形 err_ac max(abs(R_sim(:) - R_th(:))); % 零时延切片与解析式对比 C_sim abs(AF1(M1, :)).; C_sim C_sim / max(C_sim); C_th abs(sinc(T*fd)); % 需自备 sinc 函数 err_sp max(abs(C_sim - C_th));sinc在 MATLAB 里不是内置函数sin(pi*x)./(pi*x)要自己写注意x0处单独赋 1这一步写错会导致err_sp一直是 1 附近很容易被误判成「仿真错了」。两条切片的误差都应该在1e-3量级主要来自时延栅格的离散化和矩形边界的采样点截断。如果误差到了1e-1基本可以断定是索引对齐出了问题而不是公式写错。3.2 高斯脉冲二维高斯与截断的代价高斯脉冲s(t) exp(-a*t^2)的模糊函数是分离变量的二维高斯没有旁瓣|chi(tau, fd)| sqrt(pi/(2a)) * exp(-a*tau^2/2) * exp(-pi^2*fd^2/(2a))注意tau方向的指数系数是a/2而不是a比信号本身宽了sqrt(2)倍fd方向的宽度由a反比控制a越大时域越窄多普勒展宽越厉害。这就是时宽带宽积在模糊函数上的直接体现也是高斯脉冲「时频不确定性」的教学样本。但仿真里几乎不可能画出完美的无旁瓣图原因在截断。真正的高斯在±inf上才有定义代码里必然要截到有限长度。截到±2 sigma时边界处的幅度还有exp(-2) ≈ 0.135这个台阶会在频谱上产生明显的 sinc 型旁瓣画出来就是高斯主瓣周围一圈 -30 dB 左右的细纹。截到±4 sigma边界幅度降到exp(-8) ≈ 3.4e-4旁瓣压到 -70 dB 以下肉眼就干净了。for k [2 3 4] sig T/(2*k); ak 1/(2*sig^2); sk exp(-ak*t.^2); [AFk, ~, ~] ambg(sk, fs, M, Nfft); lvl 20*log10(max(abs(AFk(:)))/max(abs(AFk(M1,:)))); fprintf(截断 %d sigma边界幅度 %.2e切片旁瓣 %.1f dB\n, ... k, exp(-ak*(T/2)^2), lvl); end这段循环跑完会看到一个很直观的规律截断越狠旁瓣越高。这也是为什么做高斯脉冲的旁瓣分析时必须把截断长度写进报告否则别人复现不出你的数。3.3 三条交叉验证断言光看图不够得用解析式做断言把仿真钉死。我一般固定跑这三条断言理论值检查方式容差峰值位置(tau, fd) (0, 0)[~, imax] max(abs(AF(:)))反解行列1 个栅格零多普勒切片矩形三角形高斯exp(-a*tau^2/2)与解析式逐点求最大绝对误差5e-3能量守恒sum(sum(abs(AF).^2))dtaudfd ≈chi(0,0)^2第三条是很多人忽略的。模糊函数满足 Moyal 恒等式二维积分等于峰值模值的平方。数值上做这个积分时dtau 1/fs、dfd fs/Nfft两个步长相乘正好是1/Nfft跟采样率无关。如果积分结果偏大或偏小很多往往是AF被重复归一化了或者多普勒轴的范围取错多算或漏算了一段。dtau 1/fs; dfd fs/Nfft; E_af sum(abs(AF1(:)).^2) * dtau * dfd; E_pk abs(AF1(M1, Nfft/21))^2; fprintf(能量比 %.4f\n, E_af / E_pk); % 应接近 1注意Nfft/21这个列索引对应fftshift后的零多普勒位置因为fd(1) -fs/2。这个1忘掉的话取到的就是fd -fs/Nfft的相邻点能量比会莫名其妙地偏小几个百分点。4. LFM 信号斜脊斜率、距离-多普勒耦合与加窗代价4.1 从解析式读出脊线方程 fd K*tauLFM 的模糊函数同样有闭式解。代入s(t) exp(j*pi*K*t^2)|t| T/2在|tau| T上|chi(tau, fd)| (T - |tau|) * |sinc( (fd - K*tau) * (T - |tau|) )|跟矩形脉冲的式子只差一处sinc 的自变量从fd变成了fd - K*tau。就这一项把整张图从十字栅格扭成了一条斜线。峰值条件fd - K*tau 0也就是fd K*tau。这是一条过原点、斜率为调频率K的直线。物理含义是时延tau和多普勒fd在 LFM 上不是独立的一个tau的延迟会被匹配滤波当成K*tau的频移反过来一个多普勒频移也会被折算成fd/K的距离偏差。这就是距离-多普勒耦合耦合系数就是1/K。举个具体数B5 MHz、T20 usK 2.5e11 Hz/s。若目标多普勒为 50 kHz等效时延偏差50e3/2.5e11 0.2 us对应距离 30 米。这个量级在窄带系统里可以忽略但在大时宽带宽积的宽带系统里必须补偿。脊线的宽度沿垂直于脊的方向是1/T沿脊方向则被拉得很长这就是 LFM 对多普勒不敏感、但一旦多普勒超出脊宽就迅速掉出主瓣的原因。图上看到的那条亮线长度本质上等于[−T, T]区间内K*tau的取值范围即±KT ±B。% 用峰值轨迹拟合斜率验证脊线方程 [~, ipk] max(abs(AF3), [], 1); % 每个 tau 上的峰值多普勒索引 fd_pk fd(ipk).; mask abs(tau) 0.6*T; p polyfit(tau(mask), fd_pk(mask), 1); fprintf(拟合斜率 %.3e Hz/s理论 K %.3e Hz/s\n, p(1), K);拟合出来的斜率应该和K在 1% 以内吻合。如果符号反了说明多普勒定义里的exp(±j*2*pi*fd*t)跟你生成 LFM 时用的exp(j*pi*K*t^2)符号不匹配两处约定必须同时翻只翻一处就会得到镜像图。4.2 用 FFT 做多普勒维复杂度从 N²M 降到 N·Nfft·log第 2 章那个ambg函数其实已经用了 FFT这里把复杂度账算清楚方便判断什么时候该换算法。按定义直接算对每个tau共2M1个和每个fd共Nfft个都要做一次长度N的内积总运算量约(2M1) * Nfft * N。以N1600、Nfft16384、MN-1计约8.4e11次复数乘纯 MATLAB 循环跑要几十分钟。改成「固定tau、一次性 FFT」之后每个tau只做一次长度Nfft的 FFT总运算量约(2M1) * Nfft * log2(Nfft)同参数下约1.0e9秒级完成。差了两个数量级。还有一条路是二维 FFT把模糊函数写成s的时频分布再做二维变换复杂度还能降但实现复杂度高而且对补零和周期性假设敏感容易引入循环卷积混叠。对于N在几千以内的教学和工程验证逐tau做 FFT 已经够用没必要上二维。真正要注意的是内存。AF是Nfft x (2M1)的复矩阵Nfft16384、M1599时约 4200 万复数占 670 MB。再同时保留三个信号就是 2 GB普通机器会开始吃交换分区。我一般的做法是算完立刻归一化转成single或者只保留|AF|的 dB 矩阵。AF3 single(AF3); % 减半内存 AF3 20*log10(abs(AF3)/max(abs(AF3(:))) eps);转single之前必须先做归一化否则max在single上算出来的峰值会有相对误差动态范围 60 dB 的图上能看出 0.1 dB 的偏移做旁瓣对比时会有干扰。4.3 加窗对比Hamming 把 -13.2 dB 压到 -42 dB 的代价LFM 不加权时零多普勒切片的主瓣旁瓣比是 -13.2 dB跟矩形脉冲一样因为匹配滤波输出的包络本质上是sinc。工程上要压到 -40 dB 以下标准做法是在频域对匹配滤波器加权也就是用失配滤波换旁瓣。加权的位置在频域H(f) conj(S(f)) * W(f)输出y(tau) IFFT(|S(f)|^2 * W(f))。W(f)在信号带宽内是窗函数带宽外置零。这样做的代价是主瓣展宽-3dB 宽度从0.886/B涨到约1.30/BHamming分辨率损失约 47%。Nfft2 2^nextpow2(8*N); f (0:Nfft2-1). * (fs/Nfft2); f(f fs/2) f(f fs/2) - fs; % 搬到 [-fs/2, fs/2) S fft(s_lfm, Nfft2); W zeros(Nfft2,1); mask abs(f) B/2; W(mask) hamming(sum(m p a hrefhttps://download.csdn.net/download/zrh86/85303218 stylecolor:#ec7500;font-size:14px; 本文还有配套的精品资源点击获取 /a img altmenu-r.4af5f7ec.gif srchttps://csdnimg.cn/release/wenkucmsfe/public/img/menu-r.4af5f7ec.gif stylewidth:16px;margin-left:4px;vertical-align:text-bottom;cursor:text; /p