简介：线性调频（LFM）信号仿真与短时傅里叶变换（STFT）时频分析是雷达与通信信号处理中的基础内容。这是一个MATLAB示例，面向相关专业学生及科研人员，重点演示了LFM信号的生成、STFT时频图计算，以及通过D倍抽取提升时频分辨率的方法，适合用于理解非平稳信号分析中的分辨率权衡问题。压缩包内仅包含1个m脚本，大小约1KB，代码结构简洁，涵盖信号参数设置、窗函数选择（汉明窗等）、STFT变换、下采样抽取及结果可视化，便于直接运行和二次修改。已有881人学习，说明该示例在同类资源中具有一定参考价值。借助代码，读者可以直观看到LFM信号频率随时间线性变化在时频图上的呈现，并理解当瞬时频率变化较快导致时频分辨率不足时，D倍抽取如何改善频率分辨率，同时也需留意可能的混叠风险。对于入门雷达信号处理或开展课程仿真实验都有实际帮助。
1. LFM与STFT：为什么线性调频信号必须看时频图
雷达、声呐和通信系统里做线性调频仿真，LFM信号几乎是绕不开的波形。它表面只有一个从低频匀速扫到高频的频段，但直接用FFT看频谱，只能看到一段平坦的带宽，完全看不出频率是按什么顺序扫过去的。STFT（短时傅里叶变换）恰好补上这个缺口：把信号切成一段段加窗的短帧，逐帧做FFT，再拼成一张频率随时间变化的时频图。做LFM-STFT处理的人，大多是为了确认瞬时频率曲线是否笔直、调频斜率是否跟理论吻合、信号和干扰能不能在时频平面上分开。这张图直接决定后续测距、测速或去噪策略能否成立，所以把生成和参数调整做扎实比急着上算法更重要。
2. 把LFM信号造对：chirp参数表与仿真生成代码
2.1 线性调频信号的数学模型与三个关键参数
LFM信号的复数表达式为：
s(t) = exp(j2π(f0t + 0.5μt²))
其中f0是起始频率，μ = B/T是调频斜率，B是带宽，T是脉宽。对相位求导数，瞬时频率为fi(t) = f0 + μt。也就是说，理想时频图应该是一条斜率为μ的直线。仿真时真正要定下来的是四个数：采样率fs、起始频率f0、带宽B、脉宽T。常见做法是先定T和B，再反推μ，最后根据瞬时频率上限确定fs。
| 参数 | 含义 | 我的仿真常用取值 | 工程量级参考 |
|---|---|---|---|
| fs | 采样率 | 1000 Hz | 几百 MHz 到 GHz |
| f0 | 起始频率 | 50 Hz | 中心频率附近 |
| T | 脉宽 | 1 s | 几十 μs 量级 |
| B | 调频带宽 | 200 Hz | 几 MHz 到几十 MHz |
| μ | 调频斜率 | B/T = 200 Hz/s | 自动算 |
为什么仿真推荐用复指数而不是实数cos或sin？因为STFT的中间结果是复数谱，复数谱保留幅度和相位两个信息，后续做时频掩码、脊线提取都要用到相位。实数信号会产生正负频率镜像，时频图上出现两条对称斜线，看起来像多了一个假信号。分析阶段用复信号，到真正回放或发射再取实部，是最省事的路径。
2.2 采样率与Nyquist约束
生成LFM前必须校核一个不等式的约束条件：f0 + B < fs/2。瞬时频率从f0线性增长到f0+B，任何时候都不能越过采样率的一半，否则时频图边缘会出现回卷，斜线跑到图外又被折返进来。我在仿真前一般先用一行print语句确认上限，免得图出来之后再倒查参数。
常见误用是把信号写成sin(2π(f0+μt)t)。表面看也是“频率从f0线性增长”，实际对相位求导后瞬时频率变成f0+2μt，斜率比设计值翻了一倍。正确的相位是瞬时频率对时间的积分，也就是0.5μt²这一项，瞬时频率不能直接乘以t。这里出差错时频图很难看出来，但和理论曲线一对比就露馅。
2.3 用Python生成LFM信号的完整实现
import numpy as np fs = 1000.0 # 采样率，单位 Hz T = 1.0 # 脉宽，单位 s f0 = 50.0 # 起始频率，单位 Hz B = 200.0 # 调频带宽，单位 Hz mu = B / T # 调频斜率 = 200 Hz/s t = np.arange(0, T, 1/fs) # 时间序列 x = np.exp(1j * 2 * np.pi * (f0 * t + 0.5 * mu * t**2)) # 自检：瞬时频率必须满足 f0 + B < fs/2 fi = f0 + mu * t print(f"min freq = {fi.min():.2f} Hz, max freq = {fi.max():.2f} Hz, fs/2 = {fs/2:.0f} Hz")这段代码的核心在相位项f0t + 0.5mut**2。f0t对应起始频率的累积相位，0.5mut**2是线性调频带来的额外相位增量，这两项合在一起求导才得到f0+μt。用np.exp构造复指数信号，幅度恒为1，频率随时间线性变化，是最标准的LFM模型。
参数上，fs=1000Hz时最高瞬时频率250Hz，留一倍余量给Nyquist边界；T=1s让时频图横轴足够宽，方便肉眼对比斜率；B=200Hz让频率变化明显但不至于太陡。如果你的仿真场景是雷达脉冲，通常把T缩到微秒级、B放到MHz级，代码结构不用改动，只改这五个数值。
3. 用STFT切开LFM：窗长、重叠率与窗函数的三组实测参数
3.1 STFT的计算过程与输出尺寸
STFT做的事情可以拆成三步：把信号分帧，对每一帧加窗，再做FFT。所有帧的频谱按时间顺序堆叠起来，就是一个复数矩阵Zxx。矩阵的行对应频率，列对应时间帧，元素值表示某个时刻某个频率分量的幅度和相位。
用scipy实现只需要一行：
from scipy.signal import stft f, t_frame, Zxx = stft(x, fs=fs, window='hann', nperseg=128, noverlap=96, nfft=256)三个返回值的含义是：f是频率轴，t_frame是每个时间帧的中心时刻，Zxx是复数时频矩阵。nperseg是窗长（采样点数），noverlap是相邻两帧重叠的采样点数，nfft是FFT点数。Zxx的行数是nfft//2+1，列数由信号长度、窗长和重叠率共同决定。画图时把20*lg(|Zxx|)作为像素值，横轴取t_frame，纵轴取f，就是常见的时频瀑布图。
3.2 窗长对时频图的影响：三组参数实测对比
窗长nperseg越大，频率分辨率越高，但每个窗对应的时间跨度也越长，时间定位越模糊。这个矛盾是固定的，窗长不能同时把两个方向都调好。
| nperseg | 频率间隔 | 时间定位跨度 | 直观效果 |
|---|---|---|---|
| 64 | fs/64 ≈ 15.6 Hz | 64 ms | 斜线粗、有颗粒感 |
| 128 | fs/128 ≈ 7.8 Hz | 128 ms | 斜线平滑，边沿可读 |
| 256 | fs/256 ≈ 3.9 Hz | 256 ms | 斜线细腻，端点模糊 |
我用同一段LFM分别跑这三组参数对比过。nperseg=64时，斜线明显变宽，颜色深浅不均匀，像是锯齿状的丝带；nperseg=256时，斜线最细，但脉宽首尾各出现一段模糊过渡区，原因是窗太长，把前后的静默段也卷进来了。做LFM分析我一般从nperseg=128起步，先看整体斜率，再根据细节需求往两个方向调。需要分辨邻近两个频率分量就加长窗，需要准确定位某个时刻的频率跳变就缩短窗。
3.3 重叠率、窗函数与nfft补零的各自作用
noverlap控制帧与帧之间的重合程度。如果不重叠，相邻帧之间缺乏连续性，时频图会出现竖条纹理，斜线像被竖线切成一段一段。我一般取noverlap = nperseg * 0.75，也就是重叠加75%，这样时间方向的采样足够密集，画出来的图连续感最好。
window参数控制窗函数类型。hann窗主瓣略宽但旁瓣低，是时频分析最常用的默认选择；blackman旁瓣更低，但主瓣更宽，适合动态范围大、弱信号旁边有强信号的场景；矩形窗频谱泄漏严重，除了测试不建议用在显示上。
nfft经常被误解。nfft=256而nperseg=128时，相当于在每帧末尾补零到256点再做FFT，频率轴的间隔变小，但这只是插值，不是真正的分辨率提升。真正的频率分辨率由窗长决定，补零只让谱线看起来更密。LFM仿真里，如果要让斜线更连续，第一优先级是加长nperseg，nfft稍微大于nperseg即可，没必要设成窗长的四倍以上。
import matplotlib.pyplot as plt plt.figure(figsize=(9, 5)) plt.pcolormesh(t_frame, f, 20 * np.log10(np.abs(Zxx) + 1e-12), shading='gouraud', cmap='jet') plt.xlabel('时间 / s') plt.ylabel('频率 / Hz') plt.title('LFM信号的STFT时频图（nperseg=128, noverlap=96）') plt.colorbar(label='dB') plt.tight_layout() plt.show()画图代码里加1e-12是为了避免零幅度取对数报错。pcolormesh比imshow更适合非均匀网格，shading='gouraud'会让颜色过渡更平滑。看这张图时，重点检查三条线：左边界的斜线起点是否接近f0，右边界的终点是否接近f0+B，中间是否保持直线。如果直线弯了，回头查信号生成；如果直线断了，多半是窗长或重叠率的问题。
4. STFT处理LFM的常见坑：从黑匣子到谱图翻车的排查手册
4.1 条纹断裂和时频图掉色：窗长太短的典型症状
现象：nperseg设成32或64时，LFM斜线变成一截一截的颗粒，颜色深浅不均，斜线中间出现断续，没法准确读出斜率。
原因：窗长太短导致频率分辨率过低，一个频率bin覆盖十几甚至几十赫兹；同时窗内信号本身又在快速扫频，能量被摊到多个相邻bin上，幅度被削弱，低于colorbar显示阈值后就“掉色”了。
解决：把nperseg从64提到128或256。如果脉宽T本身很短，比如只有几十微秒，再提升nperseg会导致时间帧数太少，此时应该降低带宽B，而不是硬拉窗长。先确认频域曲线的最高频率没有越过fs/2，再做窗长调整，排查顺序不要反。
4.2 谱图回卷：采样率设错导致的“倒八字”
现象：时频图上的斜线扫到顶部边界后，没有消失，而是从底部又出现一条反向斜线，整体形成倒八字或V字形。
原因：瞬时频率f0+μt在脉宽后半段超过了fs/2，超过Nyquist的部分被折叠回低频区。这不是STFT的问题，是信号生成阶段参数不合理。
解决：首先生成信号时做硬校验，不满足约束就直接报错：
if f0 + B >= fs / 2: raise ValueError(f"瞬时频率最高 {f0+B:.1f} Hz 超过 fs/2={fs/2:.1f} Hz，必须提高采样率或减小带宽")调参数时建议先把fs提高30%以上留余量，而不是卡着Nyquist边界跑仿真。真实硬件里抗混叠滤波器有过渡带，采样率通常取最高频率的2.5倍以上，仿真可以放宽到2倍，但别贴边。
4.3 时频图只有一条亮带：动态范围太大掩盖弱信号
现象：整张图除了最亮的主斜线，周围全是黑的，肉眼看不到任何底噪，让人觉得信号特别干净。实际上噪声可能很大，只是显示范围设错了。
原因：20*lg(|Zxx|)转换后动态范围能到60dB以上，如果colorbar的vmin从最小值算起，所有比峰值低40dB的点都显示成黑色，弱信号被显示逻辑吞掉了。
解决：画图时手动限制动态范围，只显示峰值为中心的上下区间：
spec = 20 * np.log10(np.abs(Zxx) + 1e-12) plt.pcolormesh(t_frame, f, spec, shading='gouraud', cmap='jet', vmin=spec.max() - 50, vmax=spec.max())vmin取峰值减50dB，能把主斜线之外的底噪轮廓带出来。如果目标看窄带干扰或微弱的二次谐波，把动态范围压到30dB甚至20dB，视觉对比更强。这属于显示层的调整，不影响Zxx数据本身。
4.4 生成了两倍斜率的“伪LFM”
现象：时频图斜线看起来有明显斜率，但和理论瞬时频率曲线一对比，斜率是设计值的两倍，或者同时出现正负两条对称斜线。
原因：前文提到过，把相位直接写成瞬时频率乘t，即sin(2π(f0+μt)t)，求导后斜率翻倍；或者用实数信号做STFT，负频率镜像产生第二条对称斜线。
解决：统一用复指数模型，相位严格写成0.5μt²。定位这个问题最快的办法是提取脊线：找每一列幅度最大点的频率，画出来和理论fi=f0+μt对比。误差大说明信号模型写错了，误差小但图上有对称斜线说明用了实信号，改回复信号即可。
5. 时频掩码去噪：让LFM的STFT结果变得干净可用
5.1 什么时候值得上时频掩码
仿真做到中途，经常遇到两种情况：一是LFM信号被宽带噪声淹没，时频图里主斜线若隐若现；二是有窄带干扰，比如某固定频率的电台信号，在时频图上表现为一条竖直亮线，把LFM斜线切断了。这两种场景都适合用时频掩码处理。
时频掩码的思路很直接：LFM在时频平面上是一条窄斜线，能量集中，宽带噪声均匀铺满整个平面，窄带干扰只占据少数频率行。按幅度或能量设置阈值，把“不是信号”的时频点压掉，再把处理后的复数矩阵做ISTFT还原成时域信号。这套流程在语音增强里叫频谱减法或理想掩码，搬到LFM仿真的优势是不用设计带通滤波器，滤波器中心频率必须跟随时间变化，而掩码天然带时变特性。
5.2 硬掩码与软掩码：两种实现与代码
掩码分两种。硬掩码把时频点直接置0或保留，规则简单，但边界突变会在重建信号里引入“音乐噪声”，听起来像一阵一阵的啁啾声，仿真里表现为重建波形出现毛刺。软掩码给每个时频点乘一个0到1之间的权重，信号强的区域权重接近1，弱的区域权重接近0，过渡平滑，重建误差通常更小。
import numpy as np from scipy.signal import stft, istft # 构造带噪信号：原信号叠加白噪声，SNR约 5 dB rng = np.random.default_rng(0) noise = rng.normal(0, 0.3, len(x)) x_noisy = x + noise f, t_frame, Zxx = stft(x_noisy, fs=fs, window='hann', nperseg=128, noverlap=96, nfft=256) P = np.abs(Zxx) ** 2 # 时频能量 noise_floor = np.median(P) * 3 # 底噪估计，用中位数乘系数 alpha = 2.0 mask_hard = (P > noise_floor).astype(float) # 硬掩码 mask_soft = P / (P + alpha * noise_floor) # 软掩码，Wiener型权重 Z_hard = Zxx * mask_hard Z_soft = Zxx * mask_soft _, x_hard = istft(Z_hard, fs=fs, window='hann', nperseg=128, noverlap=96, nfft=256) _, x_soft = istft(Z_soft, fs=fs, window='hann', nperseg=128, noverlap=96, nfft=256)噪声底噪用np.median(P)乘以系数3，是因为LFM斜线只占时频平面的一小部分，大部分时频格子里只有噪声，中位数能代表底噪水平。系数3相当于把阈值定在底噪上方约4.8dB的位置，太大会误删信号边沿，太小则残留噪声多。alpha控制软掩码的压缩强度，alpha越大，弱信号被抑制得越狠，重建波形越干净但损伤越明显，一般取1到4之间。
5.3 掩码后的ISTFT重建约束
ISTFT能还原信号有个前提：窗函数和重叠率满足COLA条件，即各帧窗函数相加后为常数。只要生成掩码前后用的window、nperseg、noverlap保持一致，主要条件就能满足。但掩码本身破坏了每帧数据的连续性，所以重建信号会有幅度调制伪影，这是不可避免的，只能通过软掩码来缓解。
我每次做完掩码都会立刻跑一个无噪声对照：对纯净x做一模一样的掩码流程，理论上掩码全为1，重建误差应该接近浮点精度。如果这个对照都出现明显误差，问题多半在stft和istft参数不一致，比如window不同或noverlap没对应上。
另外要确认重建信号长度。istft返回的长度和原始输入不一定完全相等，因为帧数和hop数存在边界取舍。实际使用时取x_soft[:len(x)]对齐，或者用length参数强制指定输出长度。这个细节不处理，后续做相关分析时长度错位会引发莫名其妙的“时移”。
# 重建误差自检，只保留前半秒主斜线可见区域 x_soft = x_soft[:len(x)] snr = 10 * np.log10(np.sum(np.abs(x) ** 2) / np.sum(np.abs(x - x_soft) ** 2)) print(f"重建信噪比 = {snr:.2f} dB")如果SNR太低，先调alpha和noise_floor系数，再看时频掩码是否把主斜线边缘切掉。LFM的谱线周围有旁瓣，硬掩码很容易把旁瓣当成噪声干掉，导致重建信号频谱出现波纹，所以我在实弹场景里基本只用软掩码。
6. 脊线提取自检：验证STFT参数是否真的调对了
常看时频图判断“像不像直线”有一定主观性，我最后会做一步脊线提取，把视觉印象换成数字指标。方法是对STFT矩阵的每一列取幅度最大值对应的频率索引，按列拼起来就是脊线频率估计。这个频率应该和理论瞬时频率f0+μt几乎重合。
idx = np.argmax(np.abs(Zxx), axis=0) # 每列幅度最大点的行索引 f_est = f[idx] # 换成实际频率值 f_true = f0 + mu * t_frame # 理论瞬时频率 err = np.abs(f_est - f_true) print(f"平均误差 = {err.mean():.2f} Hz，最大误差 = {err.max():.2f} Hz")跑这一步时有三个判断：平均误差应该在频率分辨率的一半以内，比如nperseg=128时分辨率约7.8Hz，平均误差明显大于4Hz就要警惕；最大误差通常出现在脉宽首尾两端，因为窗在边界处覆盖不足，幅度重心偏移；如果误差持续单向偏大，回查信号模型，大概率是0.5*μt²写成了μt²。
我自己的习惯是把nperseg分别设64、128、256跑一组，记录三组的平均误差，误差收敛到分辨率极限说明参数已经调到位，再往上加窗长只会让时间维越来越糊。这套自检大概花五分钟，但能省掉后面反复看图猜参数的功夫。仿真参数这回事，图看着顺眼只能算第一步，数字说话才算落地。希望帮到你。
本文还有配套的精品资源，点击获取
