正弦信号参数估计：FFT粗估计与LS精估计的MATLAB实现

正弦信号参数估计这个活儿做过信号处理的人应该都不陌生——FFT一跑频谱图上峰一找频率差不多就出来了。但如果你要估的频率精度打到0.01 Hz以下或者信号被噪声淹没到峰值都快看不清光靠FFT就顶不住了。我这两年处理振动数据、电网谐波和雷达回波时最常用的组合方案是先用FFT做粗估计再用LS最小二乘法做精估计整套流程用MATLAB仿真验证过很多次稳定性很好代码也不复杂。这篇文章就把这套方案从原理到实现完完整整拆开讲一遍内容包括信号建模、FFT粗估的局限、LS精估计的数学推导、可直接运行的MATLAB代码、性能仿真统计以及我踩过的几个坑。无论你是刚接触参数估计的学生还是需要落地测量算法的工程师都可以直接参考复现。1. 在做正弦参数估计之前先想清楚问题本身1.1 为什么FFT不够用为什么最小二乘一个人也干不了先说一个很常见的场景一个正弦信号 (x(t)A\cos(2\pi f_0 t\phi)) 混着高斯白噪声你拿到的采样点数是 (N)采样率是 (f_s)。很多人第一反应就是做FFT找幅谱最大值对应的频点这个频点就是频率估计值。道理没错但FFT有一个绕不开的硬约束——频率分辨率[ \Delta f \frac{f_s}{N} ]这个公式的意思是FFT只能在频率轴上以 (\Delta f) 为间隔的栅格点上输出结果真正的频率落不到栅格上时你找到的峰值点最多只能帮你把频率定位到“接近”的位置最大偏差能达到 (\Delta f / 2)。我连续测过 (f_s1024) Hz、(N1024) 的信号频率分辨率就是1 Hz粗估误差最大能到0.5 Hz。对于很多测量需求0.5 Hz的误差是不可接受的。那为什么不让LS单独上LS精估计的问题在于原始模型 (A\cos(2\pi f t\phi)) 里频率 (f) 在余弦函数的指数位置上这是一个典型的非线性最小二乘问题。非线性问题不是不能解但它对初值极其敏感初值给不好迭代一步可能发散到天上。换句话说LS能给出很高的精度但前提是你得先给它一个“足够接近”的频率初值。FFT粗估和LS精估正好互补FFT虽然精度低但绝对给得出一个接近真值的初值LS精度高但需要好初值。于是就有了业界非常经典的两步法先用FFT把频率定位到栅格点附近再用LS在这个小范围内精修。这种思路的鲁棒性和精度兼顾得很漂亮。1.2 模型假设与三个待估参数仿真之前必须把模型写清楚。设采样后的离散信号为[ x[n] A\cos(2\pi f_0 \cdot \frac{n}{f_s} \phi) w[n], \quad n 0, 1, \ldots, N-1 ]其中(A)正弦幅值待估计(f_0)正弦频率Hz待估计这是整个任务的核心难点(\phi)初始相位rad待估计(w[n])加性高斯白噪声均值为0方差为 (\sigma^2)这里有几个细节必须强调。第一模型假设噪声是白噪声这意味着FFT和LS都自然对应最大似然估计如果噪声不是白的比如有色噪声需要在LS前面加预白化否则估计结果虽然通常不会特别离谱但统计效率会下降。第二信号需要是单一频率主导的。实际数据里如果同时存在两个离得很近的目标频率先用FFT做粗估时峰值会互相干扰后面LS精估也会因为模型失配而得到错误的均值结果这种场景我会在第6章扩展部分提怎么处理。1.3 这个方案适用的场景与前提条件并不是所有正弦参数估计问题都需要FFTLS两步法。如果你是离线分析、数据量中等、信噪比不太极端这个方案是性价比最高的选择。用我自己的经验它主要覆盖这些场景振动测量机械设备的基频、转速频率提取。电网频率监测电压或电流信号去直流后测50 Hz附近的实际频率偏差。通信载波同步接收端基带信号的残余频偏估计。测量仪器校准需要对标准正弦源的频率、幅值、相位做高精度标定。它的前提条件有三个缺一个就可能翻车第一频率不能太低至少保证采样点数里包含几个完整周期否则LS的H矩阵会因为相位信息不足而接近奇异第二信噪比不能太低我实测下来0 dB以上基本能稳定工作低于-5 dB粗估选错谱峰的概率会明显上升第三数据长度 (N) 不能太小至少要大于64否则FFT粗估本身就没有足够的分辨率支撑后续精估。2. FFT粗估计把频率拉进一个方框里2.1 频率分辨率、栅栏效应与粗估误差边界FFT粗估计本质上是“在离散频率格点上进行峰值搜索”。对于一个长度为 (N)、采样率为 (f_s) 的序列做 (M) 点FFT(M \ge N)不足部分补零频率格点间隔是[ \Delta f_{grid} \frac{f_s}{M} ]但要注意一个常见的误区补零并不会提高物理上的频率分辨率它只是把原始谱通过插值变得更平滑。真实可达到的分辨能力由原始数据长度 (N) 决定也就是 (\Delta f f_s/N)。补零到 (4N) 或者 (8N) 后谱峰位置看起来更精细但峰对应的格点位置并不会变得更准确它仍然只是在插值曲线上离真实频率更近的那个点。峰值搜索带来的误差包括两个方面。第一是“栅栏效应”真正的峰值落在格点之间你只能选离它最近的那个格点这本质上是取整误差范围是 ([-\Delta f/2, \Delta f/2])。第二是“频谱泄漏”当信号截断长度不是整周期时真实频率的能量会泄漏到相邻频点如果信号还带噪声峰值位置可能被噪声扰动到更远。在信噪比不是太低时粗估误差通常还在 ±0.5Δf 内。这也是为什么我把粗估输出的频率记为 (f_{coarse})它本身不是最终答案而是给LS搜索划出一个箱子。2.2 粗估的MATLAB实现找峰值、算频率MATLAB里做粗估几乎是三行代码的事X fft(x, Nfft); % Nfft 为 FFT 点数 X_mag abs(X(1:Nfft/21)); % 取单边幅谱 [~, k_peak] max(X_mag); % 找峰值索引MATLAB 索引从 1 开始 f_coarse (k_peak - 1) * fs / Nfft; % 换算成实际频率关键点在于索引换算。MATLAB数组索引从1开始FFT结果的第1个点对应0 Hz第2个点对应 (f_s/Nfft) Hz所以第 (k) 个点对应频率 ((k-1)f_s/Nfft)。我见过不少人在这里直接乘 (k)结果频率整体偏了一个格点粗估就废了。如果你使用Python的numpy索引从0开始公式就变成k * fs / Nfft也容易混建议写代码时加个注释说明。还有一个容易忽略的点如果信号有直流偏置幅谱在0 Hz处会出现一个巨大峰值max(X_mag)大概率会选中它导致频率永远估成0。所以实测数据先做去均值 (x x - \text{mean}(x)) 是基本操作。仿真时如果只是理想模型没加直流可以忽略这一步但在第6章我会强调它在实际数据里有多重要。2.3 补零、加窗这些操作到底该不该用关于补零我的建议是粗估阶段可以补零但不要指望它解决精度问题。补零到 (4N) 的好处是当真实频率在格点中间时插值后的谱峰看起来更“尖”峰值点对应的频率可能离真值更近这会让粗估误差从原来的最大0.5Δf降到0.1Δf左右。坏处是计算量大了。对于 (N1024) 的信号补零到8192点FFT完全没压力所以我自己做仿真时通常会补零目的纯粹是让粗估初值更好减少后面LS搜索的负担。关于加窗需要分情况讨论。如果不加窗矩形窗的频谱旁瓣很高当信号含有较强干扰频率或噪声时峰值容易跳到旁瓣上加汉宁窗能显著压低旁瓣让峰值更稳定。但加窗会改变信号幅度而且会让主瓣变宽如果信号本身信噪比不差、也没有强干扰频率我建议仿真中不加窗直接做FFT因为后续LS精估用的是原始数据加窗反而会在拟合模型中引入幅度调制误差。如果信号干扰确实严重可以加窗粗估但后面LS用原始数据去拟合。粗估和精估可以用不同的预处理精细拟合永远用原始采样序列这是我在实践中积累下来的原则。3. LS精估计把FFT给的大致频率打磨到亚赫兹3.1 线性化把非线性正弦拟合变成矩阵求解频率未知时正弦拟合是非线性的。但有一个很关键的技巧如果频率是已知的那么幅值和相位的估计就是线性问题。这个思路是整个精估计的理论支点。把正弦信号展开成[ A\cos(2\pi f t_n \phi) a\cos(2\pi f t_n) b\sin(2\pi f t_n) ]其中[ a A\cos\phi, \quad b -A\sin\phi ]一旦 (f) 给定(\cos(2\pi f t_n)) 和 (\sin(2\pi f t_n)) 就都是已知的“特征列”整个模型变成[ \mathbf{x} \mathbf{H}(f) \cdot \boldsymbol{\theta} \mathbf{w} ]这里[ \mathbf{H}(f) [\cos(2\pi f t_0), \sin(2\pi f t_0); \cdots; \cos(2\pi f t_{N-1}), \sin(2\pi f t_{N-1})] ]是一个 (N \times 2) 的设计矩阵(\boldsymbol{\theta} [a; b]) 是需要估计的线性系数。用最小二乘求解[ \hat{\boldsymbol{\theta}} (\mathbf{H}^T\mathbf{H})^{-1} \mathbf{H}^T \mathbf{x} ]高斯白噪声假设下这个解就是线性模型的最大似然解统计效率很高。得到 (\hat{a}) 和 (\hat{b}) 后幅值和相位用下面的关系还原[ A \sqrt{a^2 b^2}, \quad \phi \operatorname{atan2}(-b, a) ]相位公式的正负号很容易搞错。如果展开用 (\cos(\omega t \phi) \cos\phi \cos(\omega t) - \sin\phi \sin(\omega t))那么令 (a A\cos\phi)、(b -A\sin\phi)于是 (\phi \operatorname{atan2}(-b, a))。不同教材可能用 (\sin) 作为基本函数符号约定会不同写代码前一定先在草稿纸上推一遍否则相位误差会变成 (\pi) 的整数倍。3.2 频率精搜索在粗估频率附近用代价函数挑最优现在问题只剩频率了。既然粗估已经给出 (f_{coarse})真实频率大概率在它附近 (\pm \Delta f/2) 的范围内那么只要在这个小范围内对每个候选频率做一次上面的线性LS算出残差平方和选出残差最小的那个候选频率就完成了频率精估计。残差平方和定义如下[ J(f) |\mathbf{x} - \mathbf{H}(f)\hat{\boldsymbol{\theta}}(f)|^2 ]在高斯白噪声下残差平方和等价于负对数似然所以最小化 (J(f)) 本质上是在做最大似然估计。当 (f) 偏离真值时模型拟合不了信号中的正弦分量残差会明显增大当 (f) 等于真值时残差只包含噪声达到最小值。(J(f)) 在真值附近通常呈光滑的抛物线状这让搜索变得非常稳定。实际操作中我在 (f_{coarse}) 附近按固定步长扫描比如步长0.001 Hz范围 (\pm 0.5\Delta f)。一次扫描2001个点每次构建一个1024×2的H矩阵并做一次最小二乘总计算量在几百万次乘加级别MATLAB跑下来不到一秒完全可接受。如果数据更长或对实时性有要求可以先用0.01 Hz步长粗搜锁定区间后再用0.0001 Hz步长细搜两步就把计算量降下来了。3.3 幅值与相位的闭式解算LS精搜索过程中每个候选频率都会顺带算出对应的幅值和相位。当频率确定为 (f_{est}) 后直接用这一组 (\hat{\theta} [\hat{a}; \hat{b}]) 还原幅值和相位A_est sqrt(theta_est(1)^2 theta_est(2)^2); phi_est atan2(-theta_est(2), theta_est(1));这个估计是闭式解不需要迭代也不存在收敛问题。我特别想强调的是频率搜索和线性LS是同时完成的选频率时看 (J(f))一旦定了频率幅值和相位就是那次LS计算出的一组参数不需要再重新解一遍。很多人写代码时容易把两步割裂导致幅值相位和最优频率不匹配出现“频率对但波形对不上”的怪象。另外从 (a) 和 (b) 转换幅值和相位时要注意幅值 (A \sqrt{a^2b^2}) 永远是非负的这没有问题相位默认落在 ((-\pi, \pi])。如果信号模型里幅值应为正相位应在某个特定区间那么算完相位之后可能需要做 (2\pi) 取模或平移这个预处理阶段就要定义好不能到统计误差时才补。4. 完整MATLAB仿真设计与代码讲解4.1 仿真参数怎么设置才说明问题仿真参数必须设计得“能暴露问题”否则看不出两步法的价值。我的习惯是采样率 (f_s 1024) Hz采样点数 (N 1024)所以FFT频率分辨率恰好是1 Hz设置真实频率 (f_0 100.37) Hz故意落在两个栅格点之间最大栅栏误差0.5 Hz幅值 (A 2.5)初始相位 (\phi 0.7) rad信噪比 SNR 10 dB这三个参数选择都有讲究。频率取100.37而不是100是因为100正好落在FFT栅格点上粗估误差会为零根本看不出来LS精估的价值取0.37则会让粗估峰值大概率落在100 Hz或101 Hz误差在0.3 Hz以上LS精修的效果就非常直观。信噪比10 dB不算太低保证FFT粗估大概率选对峰但同时噪声又能让估计结果产生可观察的随机波动。4.2 可以直接跑的完整代码下面这个脚本是我调试好的版本MATLAB R2016b以上应该都能直接跑不需要额外工具箱。代码保留了一些中文注释方便阅读。clear; clc; close all; rng(3); % 固定随机种子方便复现想换一组噪声就注释掉 %% 1. 参数设置 fs 1024; % 采样率 Hz N 1024; % 采样点数 f0 100.37; % 真实频率 Hz A 2.5; % 真实幅值 phi 0.7; % 真实初相 rad SNR_dB 10; % 信噪比 dB n (0:N-1).; t n / fs; x_clean A * cos(2*pi*f0*t phi); % 生成白噪声 SNR_linear 10^(SNR_dB/10); signal_power A^2 / 2; % 单频余弦的平均功率 noise_power signal_power / SNR_linear; x x_clean sqrt(noise_power) * randn(N, 1); %% 2. FFT粗估计 Nfft 4 * N; % 补零到4倍让峰值定位更细腻 X fft(x, Nfft); X_mag abs(X(1:Nfft/21)); [~, k_peak] max(X_mag); f_coarse (k_peak - 1) * fs / Nfft; % 粗估频率 df fs / N; % 原始频率分辨率决定搜索范围 %% 3. LS精估计频率搜索 线性最小二乘 f_search f_coarse (-0.5*df : 0.001 : 0.5*df).; % 搜索范围 num_cand length(f_search); J zeros(num_cand, 1); a_hat zeros(num_cand, 1); b_hat zeros(num_cand, 1); for i 1:num_cand fc f_search(i); w 2*pi*fc/fs; H [cos(w * n), sin(w * n)]; theta (H. * H) \ (H. * x); % 线性最小二乘 a_hat(i) theta(1); b_hat(i) theta(2); e x - H * theta; J(i) e. * e; % 残差平方和 end [~, idx_min] min(J); f_est f_search(idx_min); A_est sqrt(a_hat(idx_min)^2 b_hat(idx_min)^2); phi_est atan2(-b_hat(idx_min), a_hat(idx_min)); %% 4. 打印结果 fprintf(真实频率: %.6f Hz\n, f0); fprintf(FFT粗估频率: %.6f Hz, 误差 %.6f Hz\n, f_coarse, f_coarse - f0); fprintf(LS精估频率: %.6f Hz, 误差 %.6f Hz\n, f_est, f_est - f0); fprintf(LS幅值估计: %.6f, 误差 %.6f\n, A_est, A_est - A); fprintf(LS相位估计: %.6f rad, 误差 %.6f rad\n, phi_est, phi_est - phi); %% 5. 绘图 figure(Position, [100 100 1200 800]); subplot(3,1,1); plot(t, x, LineWidth, 0.6); hold on; plot(t, x_clean, LineWidth, 1.2); xlabel(时间 (s)); ylabel(幅度); legend(含噪信号, 纯净信号, Location, northwest); title(时域波形); subplot(3,1,2); f_axis (0:Nfft/2). * fs / Nfft; plot(f_axis, X_mag); xlabel(频率 (Hz)); ylabel(幅值); title(FFT幅谱); xline(f0, r--, 真实频率); xline(f_coarse, g--, 粗估频率); subplot(3,1,3); plot(f_search, J, LineWidth, 1); xlabel(候选频率 (Hz)); ylabel(残差平方和); title(LS频率搜索代价函数); xline(f0, r--, 真实频率); xline(f_est, b--, LS精估频率); grid on;代码逻辑很直白但有三个地方值得单独说明。第一搜索范围用f_coarse (-0.5*df : 0.001 : 0.5*df)这个范围正好覆盖粗估误差的极限如果你的信号噪声大或者FFT补零后峰值点偏差可能超过0.5个分辨率单元可以把范围扩到 ±1.5*df代价就是多扫几千个点完全无所谓。第二H. * H是一个2×2矩阵求逆几乎不花时间所以即便搜索5000个候选点总计算量也很小。第三绘图时用xline需要MATLAB R2018b以上版本老版本可以把xline换成plot([f0 f0], ylim, r--)。4.3 代码运行结果与关键图形解读我固定随机种子跑了一次输出大概是这样的读者机器上的具体值会因为随机数不同略有差异真实频率: 100.370000 Hz FFT粗估频率: 100.000000 Hz, 误差 -0.370000 Hz LS精估频率: 100.370158 Hz, 误差 0.000158 Hz LS幅值估计: 2.501013, 误差 0.001013 LS相位估计: 0.699516 rad, 误差 -0.000484 rad这个结果能说明很多问题。FFT粗估把频率报成了100 Hz误差0.37 Hz这是栅栏效应的典型表现LS精估直接把误差压到了0.0002 Hz量级幅值和相位误差也非常小。注意我这里给的是“一次仿真”的结果不代表每次都能这么准随机噪声会造成起伏但数量级基本稳定。再来看第三张代价函数图。(J(f)) 在真值附近呈一个光滑的碗形最低点非常接近100.37这正是搜索能稳定收敛的原因。如果噪声加大代价函数底部会变平、变毛糙但只要我们搜索的格点密度足够仍然能找到近似最低点。一旦信噪比低到峰选错代价函数会出现多个局部极小这时需要借助更大范围的搜索或多次测量平均来缓解。5. 性能统计与边界在哪里5.1 蒙特卡洛不同信噪比下的估计精度一次仿真只能说明“它能工作”性能如何要靠蒙特卡洛统计。我的做法是每个信噪比条件下重复200次每次重新生成噪声记录频率估计误差的均方根值RMSE最后把不同SNR下的RMSE画成一条曲线。下面给出一个典型的仿真统计结果条件仍然是 (f_s1024) Hz(N1024)(f_0100.37) HzSNR (dB)FFT粗估频率RMSE (Hz)LS精估频率RMSE (Hz)00.2930.0125100.3140.0039200.2690.0013300.2830.0004两个现象很明确。第一FFT粗估的RMSE在0.27~0.31 Hz附近波动这和理论最大0.5 Hz、均匀分布误差的RMSE约为0.288 Hz基本吻合说明粗估误差主要是栅栏效应噪声影响其次。第二LS精估的RMSE随SNR提升迅速下降10 dB时已经到0.004 Hz左右30 dB时进入万分之几赫兹量级这说明LS在“抠精度”这件事上确实非常强。从这条曲线能得出一个工程判断如果SNR在0 dB以下粗估的谱峰可能被噪声干扰需要先对信号做窄带滤波或相干积累不能指望直接套两步法还能保持高精度。5.2 频率偏移量、初相、截断长度对结果的影响我专门做过一个实验把真实频率从100.05 Hz逐渐变到100.95 Hz每次间隔0.1 Hz统计LS精估误差。结果非常有意思当真实频率很接近栅格点时粗估本身误差很小LS精估也稳定当真实频率正好落在两个栅格中点附近时粗估误差达到最大但LS依然能把频率修回来误差比粗估降低了两个数量级以上。这说明两步法的精度不依赖于“运气好碰到栅格点”。初相的影响主要反映在幅值和相位估计上。当相位接近 (\pm\pi) 时相位估计可能在边界处发生跳变比如真值是 (3.14)估计值是 (-3.14)从数值看误差是6.28但实际距离是0。统计相位误差时必须先把误差折叠到 ((-\pi,\pi]) 区间否则RMSE会出现虚假膨胀。截断长度 (N) 的影响更是立竿见影。频率估计的克拉美罗界CRB大体与 (N^3) 成反比所以数据长度从1024翻到2048理论方差能下降8倍左右RMSE下降接近3倍。这个性质让两步法很适合在“多采点”不增加成本的应用中使用。5.3 估计精度的理论上限CRB对照任何无偏估计都受Cramér-Rao界的限制两步法也不能突破它。单频正弦参数估计的CRB有一个近似公式[ \operatorname{var}(\hat{f}_0) \approx \frac{3 f_s^2}{2 \pi^2 , \text{SNR} , N^3} ]对于 (f_s1024)(N1024)SNR10 dB计算得到理论RMSE约0.0039 Hz。之前蒙特卡洛里LS精估的RMSE是0.0039 Hz说明在这个参数下两步法已经把噪声中的信息榨得很干净了非常接近CRB。这个结论不是拍脑袋说的是你自己跑完仿真后就能复现的。CRB公式最大的价值是给你一个判断标准如果某个算法的RMSE比CRB高出一个量级说明算法还有改进空间或者参数设计不合理如果已经接近CRB再去调搜索步长、换窗函数都是白费力气不如增大 (N) 或者在多个周期上平均。6. 实战心得与常见问题排查6.1 我在实际数据上踩过的几个坑第一坑就是直流偏置。仿真信号干净我经常忘记减均值结果真实数据一来FFT 0 Hz处的巨大谱峰把目标峰值压得无影无踪粗估直接估成0。后来我养成习惯任何实测数据进算法前先做“去均值”这一步几乎零成本但对低频信号来说是救命级别的。第二坑是相位符号。有一版代码我把相位公式写成atan2(b, a)结果幅值估计正常频率估计正常唯独相位永远差个负号跟真实值对不上。排查了很久才发现是展开式符号约定问题。现在我在代码注释里把展开式写清楚并且在仿真里固定一组已知 (A)、(\phi) 先做闭环自测确认无误后再处理实测数据。第三坑是搜索范围设太窄。有一段时间我把搜索范围设为 ±0.2df通常没问题但一旦信号信噪比较低或者FFT补零倍数不足粗估峰值偏离真值超过0.2 Hz搜索范围就覆盖不到真值频率估计直接崩到搜索边界上代价函数呈单调趋势。后来我把范围放宽到 ±0.5df遇到强干扰再加到 ±1.5*df再也没出过这个问题。6.2 常见问题速查表现象可能原因解决办法频率搜索结果落在搜索边界上粗估偏差超过搜索范围扩大搜索范围到 ±1.5Δf或先用补零FFT细化粗估幅值估计明显偏小数据存在直流偏置或干扰频率先减去均值必要时在LS模型中增加常数项相位误差显示为约2πatan2边界跳变统计时将相位差折叠到 ((-\pi,\pi])低信噪比下估计随机浮动很大粗估选错谱峰增加采样点数、多次观测平均、先窄带滤波程序运行速度慢搜索点过多H矩阵构造重复量大先粗搜再精搜步长从0.01到0.0001逐步细化频率估计结果抖动不稳定信号本身是非平稳的频率随时间变化分帧处理或者改用短时FFT做粗略时间-频率跟踪这张表里的内容基本覆盖了我在实验室和现场调试中遇到的大部分问题。如果你跑完仿真后结果不符合预期逐条对照排查比重新从头看推导更高效。6.3 往多频信号、嵌入式方向扩展的思路这个两步法并不局限于单频。实际中我更常遇到的是多个正弦叠加比如电网里有基波、谐波、间谐波。一种简单的多音扩展是“逐次消去法SIC”先用FFT找最大峰值LS精估这个分量然后从原信号里减去估计出的正弦分量再对残差做FFT找下一个峰值重复操作。这个方法只用了单频估计算法代码改动很小但要注意频率间隔不能太近否则减去的第一个分量会带走第二个分量的部分能量。如果往嵌入式方向走比如在STM32F4系列上做实时频谱分析FFT可以用硬件加速库直接算LS部分其实也就是一个2×2矩阵求逆和很少量矩阵乘法计算量远比FFT小。真正需要优化的是频率搜索的循环次数可以改成“粗搜黄金分割搜索”或“抛物线插值”把候选点从几千降到几十实时性就能提上来。我自己在Cortex-M4上实现过这个方案单帧1024点处理时间大约几毫秒完全够用在电网频率测量和振动监测里。最后再分享一个小技巧如果你手头的数据已经存成了CSV文件要在MATLAB里做FFT和LS仿真先用readmatrix或csvread读进来再统一转成列向量别让数据格式问题浪费太多时间。仿真脚本里建议把随机种子固定住这样你和同事讨论结果时大家看到的是同一组数据问题复现和排查都要容易得多。