4. 基于转速脉冲的声振信号阶次跟踪实现

为什么变转速下必须做阶次跟踪

船舶推进轴系、齿轮箱和柴油机在加减速、机动航行时转速持续变化。若仍按等时间间隔采样做FFT，同一阶次分量（如齿轮啮合阶次、轴频倍频）的频率会随转速漂移，谱线被"抹宽"，出现频率模糊，幅值也被低估。阶次跟踪的核心思路是：以转速脉冲为角度基准，把等时间采样序列重采样为等角度采样序列，使固定阶次分量在角度域中频率恒定，从而在阶次谱上形成尖锐谱线。

键相信号到角度基准的转换

设键相传感器每转输出 Z 个脉冲（齿盘齿数或反光条数），第 k 个脉冲到达时刻为 t_k。相邻脉冲间转过的角度为 Δθ = 2π/Z。在 [t_k, t_{k+1}] 内假设转速线性变化（二阶精度可用三点插值），则任意时刻的瞬时转角为：

θ(t) = θ_k + (2π/Z) · (t - t_k) / (t_{k+1} - t_k),   t ∈ [t_k, t_{k+1}]

反过来，要得到等角度采样，需要求每个目标角度 θ_n = n·Δθ_s（Δθ_s 为角度采样间隔，对应每转采样点数 N_s = 2π/Δθ_s）所对应的时刻 t_n，再在原始等时间信号上插值取值。这就是"时间→角度"的重采样映射。

实现步骤与关键代码

下面给出一个自洽的最小实现：输入为等时间采样的振动信号 x、声压信号 p、键相脉冲序列 key，采样率 fs，每转脉冲数 Z，目标每转采样点数 Ns。输出为角度域重采样后的振动与声压序列，可直接做阶次谱。

import numpy as np
from scipy.interpolate import CubicSpline

def extract_pulse_times(key, fs, thresh=0.5):
    """键相脉冲上升沿时刻（秒）。key 为归一化后的键相信号。"""
    above = key > thresh
    rising = np.where(~above[:-1] & above[1:])[0] + 1
    return rising / fs

def angle_resample(x, t_pulse, Z, Ns):
    """将等时间信号 x 重采样到等角度栅格。
    t_pulse: 脉冲时刻(秒)，长度 M；Z: 每转脉冲数；Ns: 每转采样点数。
    返回 (x_angle, theta_angle)。"""
    dtheta = 2*np.pi / Z
    # 每个脉冲对应的累计转角
    theta_p = np.arange(len(t_pulse)) * dtheta
    # 目标角度栅格：覆盖脉冲区间
    theta_max = theta_p[-1]
    theta_grid = np.arange(0, theta_max, 2*np.pi/Ns)
    # 角度→时间：以脉冲时刻为节点做三次样条（角度单调，可反解）
    cs = CubicSpline(theta_p, t_pulse)
    t_grid = cs(theta_grid)
    # 时间→幅值：对原始信号做三次样条插值
    t_orig = np.arange(len(x)) / fs
    x_cs = CubicSpline(t_orig, x)
    x_angle = x_cs(t_grid)
    return x_angle, theta_grid

# ---- 示例：模拟变转速工况 ----
fs, Z, Ns = 20000, 60, 256
T = 2.0
t = np.arange(0, T, 1/fs)
# 转速从 600 rpm 线性升到 900 rpm
rpm = 600 + 150 * t / T
f_rot = rpm / 60.0
# 累计转角（积分）
theta_true = 2*np.pi * np.cumsum(f_rot) / fs
# 合成信号：1阶 + 齿轮啮合 30 阶
x = 0.8*np.sin(theta_true) + 0.5*np.sin(30*theta_true) + 0.05*np.random.randn(len(t))
# 键相脉冲：每转 Z 个，脉冲时刻为 theta_true 跨越 2πk/Z 的时刻
pulse_idx = np.searchsorted(theta_true, np.arange(0, theta_true[-1], 2*np.pi/Z))
t_pulse = t[pulse_idx]

x_ang, th = angle_resample(x, t_pulse, Z, Ns)
# 阶次谱
X = np.fft.rfft(x_ang * np.hanning(len(x_ang)))
order = np.fft.rfftfreq(len(x_ang), d=1/Ns)
peak = order[np.argmax(np.abs(X)[order > 5])]  # 忽略直流附近
print(f"主阶次峰值: {peak:.2f} 阶 (预期 30 阶)")

工程实现中的三个易错点

• 脉冲丢失/多检：齿盘缺齿或反光条污损会导致 Δt_k 异常。应先对 Δt_k 做中值滤波，剔除偏离中位数 30% 以上的脉冲，再重建角度映射。
• 插值阶次选择：线性插值在每转采样点少（如 Ns=64）时引入幅值衰减，建议 Ns≥128 时用三次样条；对冲击类信号（齿轮断齿）优先用保形插值避免过冲。
• 角度栅格与脉冲对齐：theta_grid 起点必须与某个脉冲对齐，否则阶次谱相位随批次漂移，影响声振联合的相位一致性分析。

声振联合的阶次域对齐

振动与声压共用同一组脉冲时刻 t_pulse 和同一角度栅格，因此重采样后两路信号在角度域天然对齐，可直接做阶次域互谱、相干分析或联合特征融合。这是本案例区别于单纯振动阶次跟踪的关键：角度基准是声振两路的公共时间标尺，避免了案例3中通道延迟对齐在变转速下的残余误差。