从概念到实战:随机振动分析必读)
搜索PSD前几百条基本都和“PSD文件怎么打开”“图片转PSD分层插件”有关那是平面设计圈里的Photoshop文档。但如果你做的是振动测试、结构监测、声学测量或者设备故障诊断这个缩写指向的完全是另一件事功率谱密度Power Spectral Density。这篇文章讲的是后者。我会把PSD的物理含义、数学定义、工程计算和读图技巧串成一条线。适合刚接触随机振动和动态信号分析的工程师也适合那些用了很多年采集仪器、但还没仔细想过PSD到底是怎么算出来的朋友。读完你至少能回答这几个问题为什么随机振动要用PSD而不是普通频谱g²/Hz这个陌生的单位在说什么采集仪导出的PSD曲线我应该怎么验证它可不可信1. 先解决一个最容易混淆的PSD1.1 设计师眼里的PSD与测试工程师眼里的PSD在平面设计那边PSD是Photoshop Document就是那个带图层的源文件。很多小白搜“PSD”其实是被这个问题带进来的。而在动态信号分析领域PSD是Power Spectral Density翻译过来就是功率谱密度描述一个随机信号在频率轴上能量分布的情况。这两个概念在搜索结果里经常混在一起同一个标题下可能一边是图层操作教程一边是振动试验标准。我在现场和客户对数据时遇到过几次误会我说“你把这台设备的PSD导出来看一下”对方愣了一下以为我要他保存一个图片文件。所以第一件事就是把概念隔离开这篇文章里的PSD和Photoshop没有任何关系。1.2 PSD的近亲术语ASD、加速度谱密度、自功率谱在实际工程里PSD有很多别名和相近词经常把新手绕晕。振动台随机试验里你常看到ASDAcceleration Spectral Density单位是g²/Hz或(m/s²)²/Hz它其实就是加速度信号求出来的PSD只是换了名字强调物理量是加速度。有些老规范里还叫“能谱密度”。称呼不同底层计算逻辑完全一致。还有几个容易混淆的菜单项FFT幅值谱、自功率谱、功率谱密度。很多采集软件里同时存在这几个选项导出的数据单位完全不一样。FFT幅值谱的纵轴是信号幅值自功率谱是幅值平方功率谱密度是幅值平方再除以分辨率带宽。同一个文件选错菜单出来的曲线能和别人差几个数量级。这也是我写这篇文章的原因之一PSD不是不可理解的黑盒但工程软件把概念包装得太复杂了。2. 为什么随机振动的“强度”不能直接看频谱幅值2.1 一次实测里的困惑很多工程师习惯性说的“频谱”是指对一段采集数据做FFT之后看到的幅值谱横轴频率纵轴幅值。对正弦信号这个图非常直观50Hz处冒一根尖峰高度就是正弦的幅值干净利落。但你把同一段实测振动数据导进软件加窗、FFT、看频谱很可能会发现一个问题谱线高度跟你设定的FFT点数有关。采样率不变把FFT点数从1024改成4096原来某根谱线从0.3变成了0.08而且整条曲线像一条密密麻麻的毛刺带不会收敛到某个固定形状。这不是仪器坏了也不是信号调理出问题。随机振动本质上是非周期的连续信号它的能量分布在连续的频率区间上而不是集中在少数几条离散谱线里。你拿FFT逐点去量它量到的其实是“某一小段频率宽度内的总能量”频率宽度一变读数自然变。2.2 用“水桶接瀑布”理解分辨率带宽我当时给新同事打过一个比方想象你在瀑布下面摆一排水桶去接水每个桶代表FFT里的一根谱线桶的宽度就是分辨率带宽RBW。桶越宽每个桶接到的水越多但瀑布总的“流量密度”没变。桶越窄每个桶接到水越少但如果你算单位宽度内的水流量结果基本是稳定的。我们关心的物理量应该是“单位宽度内的流量”对应到信号里就是“单位频率内的功率”也就是功率谱密度。这就是分辨率带宽RBW这个隐藏变量的作用。FFT的谱线间隔是Δf fs / Nfs是采样率N是FFT点数。对正弦信号能量集中在一个点上无论Δf怎么变幅值基本不变对随机信号Δf变窄落在单根谱线里的随机能量变少幅值就跟着变低。所以直接拿FFT幅值谱判断随机振动大小结论不可复现。PSD解决的就是这个可复现性问题。它把每个频点上的功率再除以Δf得到“单位频率宽度内的功率密度”。你把Δf从1Hz改成0.5Hz单根谱线代表的面积变小了但除以更窄的带宽后密度值基本不变。这就是为什么随机振动标准里给出的试验条件都是PSD曲线而不是普通的频谱幅值。2.3 FFT频谱和PSD的使用边界维度普通FFT幅值谱PSD功率谱密度纵轴含义频率成分的幅值单位频率带宽内的功率密度典型单位V、g、m/s²、PaV²/Hz、g²/Hz、Pa²/Hz对离散正弦幅值准确峰值不随RBW变化显示为窄尖峰但工程上读正弦幅值一般不直接用PSD对随机信号幅值随机波动随RBW减小而变小曲线平滑基本不随RBW变化适用场景正弦成分检测、谐波分析随机振动、宽带噪声、能量分布评估FFT幅值谱并不是一无是处查谐波、找离散干扰频率它更直观。但评估随机振动能量必须用PSD。这两个工具是配合关系不是替代关系。3. PSD的数学定义和它背后的物理直觉3.1 维纳-辛钦定理从自相关到功率谱数学上功率谱密度有几种等价定义。对平稳随机信号x(t)先定义自相关函数Rxx(τ) E[x(t)x(tτ)]它描述信号和自身延迟τ之后的相似度。维纳-辛钦定理说功率谱密度就是自相关函数的傅里叶变换Sxx(f) ∫ Rxx(τ) e^(-j2πfτ) dτ这个式子很简洁但第一次看的人通常会问为什么绕这么大一圈直接对信号做FFT不行吗原因在于随机信号的样本不确定性。每一个具体的时域样本做FFT结果都带着随机起伏不是一个稳定的统计量。而自相关函数把随机性做了期望运算再用傅里叶变换拆到频域得到的Sxx(f)就是统计意义上的功率密度。另一种更贴近工程的定义是Sxx(f) lim(T→∞) E[ |X_T(f)|² / T ]其中X_T(f)是长度为T的信号段做傅里叶变换的结果。除以T是为了把总能量归一化成平均功率取期望是为了抑制随机起伏。这个定义后来在Welch法里被直接实现工程软件里跑的其实就是它。3.2 双边谱与单边谱那个2倍关系不能忘理论推导出来的通常是双边PSD频率从负无穷到正无穷都有定义而且功率是对半分到负频率和正频率的。实际工程只关心正频率所以把正频率部分乘以2得到单边PSD这是仪器和软件默认给出的结果。如果你在软件里把单双边选项选错了后面算总均方根值会直接差√2倍因为功率差2倍均方根是功率开根号。遇到过不少朋友拿着数据来问我为什么算出来和仪器对不上最后发现就是这里选错了。所以拿到陌生软件时第一件事是确认它输出的是单边还是双边PSD。3.3 单位里藏着的信息g²/Hz到底是什么PSD的单位是“信号物理量的平方除以频率”。加速度信号求出的PSD单位是g²/Hz电压信号是V²/Hz声压信号是Pa²/Hz。以加速度为例如果一条PSD曲线上某个频点的值是0.01 g²/Hz意思是在这个频点附近的1Hz带宽内振动的平均“强度”是0.01 g²。注意g在这里是重力加速度单位不是克。1g等于9.80665 m/s²所以g²/Hz换算成(m/s²)²/Hz时要乘以9.80665²大概是96.2。这是一个高频翻车点。有人把g²/Hz换算成(m/s²)²/Hz时只乘了9.81结果全线偏小一个数量级。关键就是g已经带了平方换算系数也要平方。3.4 把PSD理解成“浓度”比理解成“功率”更顺手对很多人来说“功率谱密度”里的“功率”两个字反而碍事因为测加速度时根本没有功率这个概念。我更喜欢把PSD理解为“频谱浓度”它告诉你振动能量在频率轴上怎么分布哪一段浓度高哪一段浓度低。浓度和总量是两回事。一杯浓盐水体积很小总的盐量可能不如一大桶淡盐水。对应到PSD窄峰虽然“浓度”谱密度值高但占的频率宽度窄总量未必大宽缓的凸起虽然峰值低但占据几十甚至上百赫兹总量可能远超那个尖峰。这个直觉在后面读曲线时会非常有帮助。4. 用Welch法算PSD时那些容易被忽略的工程细节4.1 Welch法四步走实际工程里最常用的PSD估计方法叫Welch法又叫改进周期图法。它把长信号分段、加窗、做FFT、平均有效降低随机误差。步骤是把长度为N的数据分成M段每段长度L相邻两段重叠D重叠率通常取50%或75%。每段乘以窗函数比如汉宁窗。对每段做L点FFT取幅值平方按窗函数能量修正得到一段“周期图”。对所有段的周期图做平均得到最终PSD。为什么必须平均因为随机信号单次FFT的幅值有极大的统计起伏只看一段谱线会像锯齿一样乱跳。平均之后统计方差大概反比于参与平均的段数。段数越多曲线越平滑结果越稳定。Welch法的代价是频率分辨率变差。分段越短RBW fs / L就越大相邻两个频率成分越难分辨。这是个绕不开的权衡想要高分辨率就要用长段但段数变少曲线变毛想要平滑稳定就要用短段但会把邻近的窄带特征糊在一起。4.2 窗函数怎么选矩形、汉宁、平顶加窗是为了减轻频谱泄漏但每种窗都有自己的取舍。做随机振动PSD时我几乎只用汉宁窗默认重叠率75%兼顾分辨率和曲线平滑度。窗主瓣宽度旁瓣抑制适用场景矩形窗最窄很差旁瓣高瞬态信号、校准信号检测优先保分辨率汉宁窗较宽较好随机振动和一般连续信号PSD的首选平顶窗最宽很好正弦幅值精确测量比如传感器标定平顶窗在随机PSD里不常用。它的幅度平坦特性适合精确读正弦幅值但主瓣太宽会把相邻频率成分抹在一起而且噪声带宽更大。如果拿平顶窗去做随机振动PSD得到的曲线会明显变“钝”共振峰被展宽边缘细节丢失。4.3 加窗之后的修正系数不能忘窗函数不只是改变泄漏它还会改变信号的总能量。计算时不能只拿|FFT|²就直接用必须除以窗函数的能量修正系数S2 Σ w[n]²否则结果偏小。标准做法是P_segment[k] |X_w[k]|² / (fs · S2)其中X_w[k]是加窗后该段数据的FFT复数结果fs是采样率。得到每段周期图后对所有段做平均再对非直流、非奈奎斯特频率乘以2得到单边PSD。这个修正系数在自写算法时特别容易漏。很多人自己写了个FFT功率谱结果比商业软件小一截多半就是漏了S2或忘了除以fs。我建议能用现成的库函数就不要自己手搓但自己要能验证库函数的输出是否合理否则出了问题根本不知道去哪找。4.4 分辨率带宽和统计平滑度之间的权衡频谱分析仪上有个参数叫分辨率带宽RBW实时频谱分析里它由FFT点数决定RBW fs / L。RBW越宽落入单根谱线的噪声越多平均后的PSD曲线越平滑但这不是“分辨率高”的表现反而是牺牲了对相邻窄带特征的辨识力。RBW越窄曲线毛刺越多统计方差越大总体趋势更难看清。工程上要反过来设计先确定需要分辨的最小频率间隔再决定段长。比如分析轴承故障特征需要分辨转速频率附近几赫兹的边带RBW就要取1Hz甚至更低若只是评估总体随机振动量级用4Hz、8Hz的RBW即可。这个取舍直接决定你最后能不能看到故障特征频率。5. 读PSD曲线的实战经验单位、总均方根和误判5.1 把总均方根算出来不要被曲线形状骗了设备厂家给的振动指标通常是总均方根值单位是g或mm/s。PSD和总RMS的关系很简单RMS_total sqrt( Σ ( PSD(f_i) · Δf ) )操作上就是先把每个频率点上的密度值乘以频率间隔Δf得到该频段内的功率分量全部累加后再开根号。这个计算在Excel里就能做也可以用代码跑。有一次客户发来一份采集仪导出的PSD对方说“峰值看起来不大设备应该没事”。我算了下总均方根0.85g已经快接近设备报警线。原因是能量主要聚在一个很宽的低频凸起里单看谱峰高度根本看不出总量。所以读PSD的第一件事不是看峰有多高而是把总面积算出来。5.2 常见的单位换算翻车现场工程标准里PSD常用g²/Hz偶尔也用(m/s²)²/Hz还有用dB表示的。三者之间的换算很容易出错。先看g²/Hz和(m/s²)²/Hz1g²/Hz ≈ 96.2 (m/s²)²/Hz。要注意是乘96.2不是乘9.81。再看dB形式。很多振动标准里写“dB re 1 g²/Hz”意思是以1 g²/Hz为基准的分贝数。换算公式是10·lg(数值)不是20·lg。20·lg是电压或幅值用的功率密度是功率量纲用10·lg。比如0.04 g²/Hz换算成dB就是10·lg(0.04) ≈ -14dB。如果把lg前面的系数用错结果会差一倍这在校对报告时非常致命。5.3 峰值高不等于能量大凸起宽不等于能量小再强调一遍浓度和总量的关系。窄尖峰看起来吓人但可能只有0.2Hz宽对总RMS的贡献其实有限。反过来一个从50Hz延伸到200Hz、峰值只有0.005 g²/Hz的宽缓凸起面积算下来可能是总RMS的大头。这种宽峰通常是宽带随机激励下的结构共振或者是某种涡激振动单看峰值很容易低估。另外如果PSD上出现异常尖锐的峰且频率刚好对应转频或啮合频率先别急着下结论说振动超标。要确认它是真实机械振动还是电气干扰、传感器谐振或者FFT处理的泄漏伪影。我一般会结合时域波形、相干函数和倒频谱一起判断绝不单凭一幅PSD下结论。6. 工程应用里PSD是怎么被直接使用的6.1 随机振动试验的条件传递和控制闭环做环境可靠性试验时随机振动条件通常用若干段PSD折线描述。比如20Hz处0.01 g²/Hz80Hz处0.04 g²/Hz然后以多少dB/oct的斜率滚降。这一组数据就是振动控制仪要复现的目标谱。控制仪的工作逻辑是采集振动台面的加速度信号在线计算PSD与目标谱比较计算误差再调整驱动信号形成闭环迭代。试验前后通过对比台面和产品上的实测PSD可以判断产品有没有被结构放大或者安装方式是否引入了额外共振。在这个流程里PSD不是事后分析的附加状态它就是控制回路里的核心反馈量。6.2 环境振动与人体舒适度评估评估办公楼、天桥这类结构物的环境振动时经常把加速度PSD按ISO 2631的频率加权函数做加权再换算成加权加速度均方根。这个方法背后的逻辑是人对不同频率的振动敏感度不一样4到8Hz的垂向振动最敏感低频段和高频段相对迟钝。如果你手里有原始PSD频域加权可以直接实现在每个频率点乘以对应的权重系数再按总RMS公式计算。这样就不用把信号分帧、滤波、再加权那么费劲了。实际项目里很多环境振动评价报告输出的就是加权后的RMS值但底层的频域代价来自PSD。6.3 故障诊断和声学里的PSD使用边界轴承早期故障的特征频率峰值微弱而且故障冲击能量分散在宽频带上直接对原始信号做FFT往往看不到明显尖峰。更常见的做法是把信号在时域做带通滤波、包络解调然后再对包络信号求PSD或FFT在包络谱里找故障特征频率。这里PSD的价值不在于直接暴露故障而在于提供了一种统计稳定地度量微弱随机信号的方法。声学里的应用更直观用声压PSD评估噪声结合1/3倍频程分析判断噪声能量集中在哪个频带进而定位噪声源。比如风机噪声的中低频凸起往往对应气动噪声高频尖峰可能对应叶片通过频率。所有场景的共同点是面对随机、宽带信号PSD是可比、可复现的度量面对主导性的离散正弦成分窄带FFT或跟踪滤波器反而更直观。选哪个工具取决于你面对的信号性质。7. 用一段Python代码把PSD算出来7.1 模拟一段加速度信号直接贴采集仪的真实数据不方便我先生成一段仿真信号模拟一个受随机激励的谐振系统再叠加一个50Hz工频干扰物理量按照重力加速度g来算import numpy as np from scipy import signal fs 5120.0 # 采样率 5120 Hz N 5120 * 120 # 120 秒数据 rng np.random.default_rng(42) t np.arange(N) / fs # 随机激励通过一个20~200Hz的带通共振系统 b, a signal.butter(2, [20 / (fs / 2), 200 / (fs / 2)], btypeband) x signal.lfilter(b, a, rng.standard_normal(N)) # 叠加一个50Hz小正弦模拟电气干扰 x 0.05 * np.sin(2 * np.pi * 50.0 * t)这个仿真信号里随机成分占主导在20到200Hz之间有明显的能量集中另外在50Hz处有一个小的离散尖峰。7.2 用Welch法求PSD并验证直接调用scipy的welch函数这是经过大量用户验证的成熟实现f, psd signal.welch( x, fsfs, nperseg4096, noverlap2048, windowhann, scalingdensity )设窗口长度4096点时RBW fs / nperseg 1.25Hz。这个分辨率足够把20到200Hz的共振峰和50Hz的离散尖峰分离开但如果要分析更精细的边带比如间隔不到1Hz的两个相邻峰就需要把nperseg加大到16384甚至更大。代价是平均次数减少曲线会变得更毛糙。7.3 计算结果的自查清单我自己处理PSD数据时每次都会做三个自查能拦住大部分数据质量事故。第一个是纯正弦校验。仿真信号里50Hz正弦的幅值是0.05g它的功率贡献约为0.05²/2 0.00125 g²。这些功率落在1.25Hz的带宽里PSD峰值应该接近0.001 g²/Hz。如果你的计算结果在这个数量级附近说明单位、单双边、窗修正基本没搞错。第二个是总均方根比对。用代码把PSD面积开根号df f[1] - f[0] total_rms np.sqrt(np.sum(psd * df)) print(total_rms) print(np.std(x))把求出的总RMS和时域信号直接算标准差的结果对比应该在同一量级。如果差异超过20%先检查是不是单边双边选错或者单位换算出了问题。第三个是参数留痕。采样率、FFT点数、窗函数、重叠率、平均次数这些参数必须写进实验记录。PSD曲线如果不记录这些参数等于没有报告因为别人无法复现你的结果你也无法在几个月后重新对比评估。我见过太多只留一条曲线、不留参数的报告最后成了谁也说不清的悬案。这三个自查做完手里的PSD数据才真正可信。