1.1 功率谱密度（PSD）简介
1.1.1 功率谱密度（PSD）
为何功率谱密度（PSD）重要
一个随机过程或随机波形本质上是无规则、非确定且非周期的。观察者无法通过测量一个时刻的数据来推测前后时刻对应的数据。随机过程或随机波形可以从统计观点来分析与理解，但是对任一点的数据，在其发生之前我们一无所知。
图 1.1 加速度-时间随机时程图（“时程”指“时间历程”，后同）
随机波形在振动测试中是一种常用的测试模式，因为它们直接与现实世界的振动相关。世界上发生在我们周遭的事件与行为都是随机的，对于即将在现实世界中使用的产品来说，贴近实际的测试是重要的。
分析随机振动最常用的工具是功率谱密度（Power Spectral Density，缩写为 PSD）。在功率谱密度图中，隐藏在时程图中的共振与谐波能够清晰地显现。在实践中，检查与分析一个随机波形的第一步通常是生成它的 PSD 图。
然而，生成一个有用的 PSD 图并非易如反掌。测试工程师和技术人员需要调节各种参数来定制 PSD 的计算以适应特定的测试情况。为了使得这些参数设置是有效的，测试工程师和技术人员需要对 PSD 的内涵有清晰的掌握。
接下来的两节课
接下来的两节课对 PSD 计算的多个步骤给出了文字表述，并给出了示例图片。尽管每一步都是一个数学推导的过程，但这些表述在阐释关键概念的同时回避了公式。如果需要了解关于 PSD 在数学上严格的表述以及相关的计算，请看1.2 节: 记录数据的转换。
1.1.2 名称中三个词的含义
“功率谱密度（Power Spectral Density）”中的每个词都代表一种基本的分析组件。
“功率（Power）”
指的是 PSD 的大小是被分析信号的均方值这一事实。它不涉及物理量“功率”(如以瓦特或马力为单位)。但是，由于功率与某个量 (如电路中的电流或电压) 的均方值成正比，任何量的均方值都被称为该量的功率。
“谱的（Spectral）”
指的是 PSD 是频率的函数这一事实。PSD 表示信号在频率谱[1]内的分布，就像彩虹表示光在波长（颜色）谱内的分布一样。
“密度（Density）”
指的是 PSD 的大小被归一化为单位赫兹带宽这一事实。例如，对于单位为 G 的加速度信号，对应的 PSD 的单位为 G[2]/Hz。
由于 PSD 这个名字不包含被测量的量，所以“功率（power）”一词有时会被被测量的量的名称所代替。例如，加速度信号的 PSD 有时被称为加速度谱密度。
1.1.3 生成一个 PSD
Vibration Research 软件采用 Welch 法进行 PSD 估计。步骤如下：
♦ 这一过程开始于正态的时域的输入数据（一个时程文件）
♦ 这数据被分成时间上等长的帧。每一帧都采用快速傅里叶变换（FFT）转换到频域上。
♦ 通过将每个频率点的量值进行平方将复杂的频域数据转换为功率。对每一帧平方后的量值
（功率值）累加后取平均值。
♦ 除以采样率以归一化到单位 Hz。
一些其他要点特别是加窗和重叠会影响最后的总体结果。
让我们看看每一步的细节。
高斯随机数据
下图给出了 5 秒钟的高斯随机数据——一个给出加速度-时间2关系的时程文件。除了能看出最大加速度大概在-30G 左右以外，从图上我们很难提取出有意义的信息。
[1] “谱”的概念最早源自光学研究，指连续的某个量值的范围。“频谱”指频域上某段连续范围。后文“信号的频谱”这种表述指信号在频谱上的幅值分布曲线。——译注
[2] 原文这里是频率，根据上下文以及配图，这里最好还是使用“时间”——译注
图 1.2 时长 5 秒的时程图
转换到频域
为了得到关于这个数据有用的信息，必须在频域里进行观察。这里用到两种计算，FFT 和 PSD。
快速傅里叶变换（FFT）将数据转换为加速度-频率图。在观察实时数据或者处理时程文件时，常用
FFT 图来监控频谱并关注频谱的变化。在生成 FFT 图时不需要窗函数、平均函数或者归一化函数。
然而，要看到频域上的能量分布，PSD 的计算是必要的。
功率谱密度（PSD）
- 计算始于将时程文件分成等时间长度的帧。
用谱线和采样率来确定每帧的宽度。每个谱线两个样本。下图时间历程的采样频率为 8192Hz，有 4096 条谱线，所得帧宽为 1 秒。如果谱线数目变为 1024 条，帧宽将会是 0.25 秒。在选择采样率的时候，这是需要牢记的重要的一步。采用 2 的幂（2n）的采样率常能使所得 PSD 的线具有便于处理的间距。
图 1.3 时程图被分成若干帧
- 在施加窗函数后，对每一帧做 FFT。
FFT 假定数据是无穷序列，这意味着每帧被程序解释为是首尾相连的。对于随机数据而言，往往不会有这种情况，因此需要施加一个窗函数。没有窗函数，每帧的首尾两点可能是不同的，从而在这两点间会得到一个瞬态峰。两个样本间的瞬态峰会在 FFT 中表现为高频能量。设想一个终端峰值冲击脉冲，从峰值到零加速度的急剧转变产生的高频能量比像半正弦这样的平滑脉冲要大得多。对于 FFT 而言，首尾点间的急剧过渡是不连续的。这种不连续性反映在 FFT 计算中，称为频谱泄漏。窗函数可以消除不连续性的凸显，减少频谱泄漏。在理想情况下，真实数据在每一帧中都有相同的起点和终点。由于情况并非如此，我们必须通过使用窗函数来将这种不一致带来的影响降到最小。
图 1.4 加窗后数据的帧
每种窗函数都有其自己的特性以适应具体的应用场景。一个窗函数由两个关键组分来估计 ——副瓣和主瓣。在上图中采用了汉宁（Hanning）窗。它有一个很高很宽的主瓣和一个几乎为零的很低的副瓣。这意味着在首尾点间几乎没有间断，从而得到非常精确的频率测量。灰色波形为原始波形，橙色波形为加窗后的数据。
汉宁窗或布莱克曼（Blackman）窗函数是振动测试和分析中最常用的窗函数，因为他们具有良好的频率分辨率和最小的不连续性，所以也就有最小的泄漏。其他窗函数可能适用于其他应用场景。想了解所有常用的不常用的窗函数的描述及其性质，请参看 VR 大学课程：用于信号处理的窗函数。
每帧加窗后的数据用来计算 FFT——一个将信号从时域向频域转换的数学函数。这种线性变换提供了观察时程波形的频率内容的能力。
图 1.5 计算每一帧的 FFT
FFT 在信号分析和处理中有许多应用场景。最重要的是，它显示了波形中的频率和比例。这可以用来确定在一段时间内什么频率被激发、一个加窗帧内每个频率的峰值加速度、峰值的分布、谐波含量等。
有一些加权因子可以应用于 FFT。如果这个过程的目标只是生成一个 FFT，那么可以应用一个加权因子来确保时域数据的 1Gpk 等于 FFT 中的 1Gpk。当目标是生成 PSD 时，对窗函数进行归一化以保持输入功率不变。
- 将每帧 FFT 平方后取平均。
术语“功率”在电气工程中是一个常见的术语，常用于指代任意值的幅值的平方。
图 1.6 对各帧的 FFT 平方与取平均
PSD 显示了一段时间内单个频率的平均能量。它最初会方差很大或很粗糙[1]；随着更多含有相似数据的帧被包含在平均值中，总体方差会减小，精度会提高，PSD 看起来会更加平滑。PSD 中包含的总时间与平均参数——自由度 (DOF) 有关。自由度越高，平均在一起的数据的帧数就越多。更多的自由度也需要更长的时间来获取和计算。
- 将计算结果归一化到单位赫兹。
这个过程的最后一步是取平方后平均的 FFT(振幅 2) 除以采样率。这使一切都归一到单一赫兹，并得到了一个功率谱密度。对于加速度而言，其 PSD 的单位是 G2/Hz。
图 1.7 计算得到的 PSD
使用 PSD 使得被测产品的反应清晰明了。随着 PSD 中数据帧数的增加，方差会不断减小，
PSD 会变得更加平滑。
要更深入地了解方差和平滑 PSD 的方法，请观看关于瞬时自由度的网络研讨会，这是 Vibration Research 的专利成果，可以快速有效地减少方差并创建非常平滑的 PSD。该方法是在短时间内显示平滑 PSD 轨迹的唯一在数学上证明了的方法。
在 PSD 创建期间，常用到另一种非必需的技术。
重叠用于在 PSD 中包含更多的原始数据，并在一段时间内生成更多的 DOF。在 0% 重叠的情况下，每一帧数据都是完全分离的。这意味着由于施加了窗函数，在每一帧中都有一些数据没有被考虑进去。此外，对于 PSD 中包含的每一帧数据，计算总平均值的两个自由度。
这意味着，如果我采用一个汉宁窗函数创建一个零重叠, 120 DOF, 8192 Hz 采样率，4096 行分辨率
(结果帧宽为 1 秒) 的 PSD，我需要对 60 秒的数据做平均来实现 120 个自由度。
50% 的重叠是指在每帧的起始点的间隔仅为 0.5 秒。每帧长度是一秒。当帧重叠发生的时候，每个 FFT 就不再有两自由度。50% 的重叠会导致每个 FFT 约 1.85 个自由度；75% 的重叠会导致每个 FFT 约 1.2 个自由度。因此，对于 50% 重叠的 120 自由度、8192Hz、4096 分辨率谱线数的汉宁窗 PSD，我能够得到期望的 64.8 帧的 PSD，但是它仅需 32.4 秒就能获取到。
图 1.8 重叠
在 0% 重叠的初始案例中，有五帧数据，导致了 5 个 FFT。如图1.8所示，具有 50% 重叠的同样一组数据得到 9 个 FFT。
PSD 计算中考虑的最后一个参数是“分辨率谱线数”。这一参数和采样速率一样决定着 PSD 中两个分析点的间隔。更多的线数会导致更精确的 PSD，但是需要更多的样本数目来进行合适的计算。
有许多测试标准要求在谐波中包含一定数量的谱线，以正确显示峰值。如果 PSD 中使用的谱线数太少，其结果类似于对波形欠采样时的情形; 分析点之间的距离太大，没有适当地考虑到它们之间的间距。一个谐波的合理采样至少需要 3+ 条谱线。
1.2 记录数据的转换
1.2.1 PSD 是什么
在振动分析中，PSD 表示信号的功率谱密度。每个单词都被选择用来表示 PSD 的一个基本组件。
“功率（Power）”指的是 PSD 的大小是被分析信号的均方值这一事实。它不涉及物理量“功率”(如以瓦特或马力为单位)。但是，由于功率与某个量 (如电路中的电流或电压) 的均方值成正比，任何量的均方值都被称为该量的功率。
“谱的（Spectral）”指的是 PSD 是频率的函数这一事实。PSD 表示信号在频率谱内的分布，就像彩虹表示光在波长（颜色）谱内的分布一样。
“密度（Density）”指的是 PSD 的大小被归一化为单位赫兹带宽这一事实。例如，对于单位为 G 的加速度信号，对应的 PSD 的单位为 G2/Hz。
由于 PSD 这个名字不包含被测量的量，所以“功率（power）”一词有时会被被测量的量的名称所代替。例如，加速度信号的 PSD 有时被称为加速度谱密度。
为何是功率？
均方值（功率）是信号强度的一种方便的度量。这在图1.9中得到说明，图中显示了加速度计测得的一辆汽车的底盘的振动时程。信号的平均幅值不能由均值给出，因为均值接近零。于是转而将信号平方（得到正值）以后再取平均。为了得到一个一次的值（图中情况下为以 G 为单位），对结果又取了方根值，从而得到 RMS（均方根）值。
图 1.9 汽车底盘的振动（均方值为 0.0053 G2，RMS 值为 0.073 G）[2]
单位正弦波的均方值为 0.5，RMS 值为 0.707。当两条正弦波叠加后，均方值为 1，RMS 值为 1。
图 1.10
从数学上讲，这基于两个独立变量 A 和 B 的一般结果。和的平方为：(A + B)2 = A2 + 2AB + B2。如果两个变量相互独立并且平均值为零，则 2AB 的平均值就为零，正如图1.11所给出的一样。
图 1.11
应该注意的是，在某些情况下，PSD 是由计算值的平方根表示的，所以加速度 PSD 的单位可能是G/(Hz)1/2。确定 PSD 的单位时必须谨慎。
为何是谱的？
当处理含有谐振的系统时，信号的频率分布是很有用的信息。这在图1.12中得到说明，其中一个悬臂梁在基础位置由一个宽带信号（含有较宽的频率分布范围）驱动，一个加速度传感器在端部测量。从信号时程中很难确定梁的谐振频率，而端部振动频谱的峰值点却清晰地显示了谐振频率。
图 1.12
为何是密度？
信号的频率分布的幅值取决于分布中频带的数目。这在图1.13中得到说明，其中采用了三种不同的带宽计算汽车振动信号的频谱。谱的幅值的平方与频率带宽成比例。为了克服这种变动，PSD 将幅值的平方除以频带宽，从而得到一个跟所取频带宽无关的一致的值。
图 1.13
需要注意的是，在某些情况下，与带宽相关的频谱和 PSD 是混用的。必须注意确定是否使用带宽归一化了幅值——这是 PSD 计算中的一部分。
1.2.2 PSD 如何计算
计算信号的 PSD 的基本步骤是确定频谱。这是通过傅里叶级数（以 1807 年发明这种方法 [1] 的数学家约瑟夫·傅里叶命名）来实现的。它是基于这样一个原理，即任何长度为 T 的信号 x 都可以由一系列的正弦和余弦波的和构成，每个波的周期数都是 T 的整数倍。数学上表述为：
(1.1)
其中 f1 = 1/T，fi = if1，Ci 和 Si 是需要确定的傅里叶系数。（x¯ 是 x 的平均值），这点在图1.14中得到形象的说明。
图 1.14
在实践中，信号 x 被数字化为采样点间隔为 ∆t 的 N 个数的数列 xn，n 取 1 到 N，如图1.15所示。
图 1.15
每项傅里叶系数的值由式1.1的两边同乘以相应的正弦或余弦波并计算对于 N = T/∆t 的均值得到。
例如，C37 由下式求出：
(1.2)
其中求和计算了均值，并且 tn = n∆t。
这一操作是可行的，因为不同的正余弦函数相互独立。即，任何两个不同正余弦函数乘积的均值为零，如图1.11。式1.3右侧除了含有 C37 的项之外所有项的均值都为零。从而：
(1.3)
其中用到了余弦函数的均方值为 1/2 这一事实（如图1.10所示）。一般地，傅里叶系数由下式给出：
(1.4)
历史地看，式1.4中求和的计算在以前是很费时的，直到 1965 年 Cooley 和 Tukey 为此提出一个有效的算法 FFT（Fast Fourier Transform，快速傅里叶变换）[2]。他们的算法采用 N 个数字采样序列，其中 N为 2 的幂（如 2,4,8,16,32 等）。FFT 可以解析的最高的频率是 fq = 1/(2∆t)Hz，称为奈奎斯特（Nyquist）频率，故而共有 N/2 个频率值。信号 x 必须在被数字化以前采用低通滤波以保证高于 fq 的频率都被去除。
FFT 计算了傅里叶系数 Ci 和 Si。均方幅值由 Ci2 + Si2)/2 给出。由于频率分辨率为 ∆f = 1/T =
1/(N∆t)，通过用均方幅值除以 ∆f，PSD 归一化到 1Hz 带宽：
(1.5)
作为 PSD 创建过程的一部分，必须处理许多 FFT 的计算细节。这些 FFT 的细节包括数字时间采样、混叠、加窗、频率分辨率和平均。
参考文献
- Fourier, J.,“Mémoire sur la propagation de la chaleur dans les corps solides,”Présenté le 21 décembre 1807 à l’Institut national –Nouveau Bulletin des sciences par la Société philomatique de Paris I (6), 112–116, 1808.
- Cooley, J. W., Tukey, J. W., “An algorithm for the machine calculation of complex Fourier series,” Comput. 19, 297–301, 1965.
欢迎关注VR公众号