汽车RPM信号FFT分析：自定义轴单位与商用级参数配置实现

问题说明 #

我尝试对汽车RPM信号执行FFT转换并绘图，但现有代码无法实现商用工具级别的频谱分析需求：需要将Y轴单位设为[n/rpm]，同时支持重采样2048Hz、频谱大小8192、汉宁窗、50%重叠、峰值幅度缩放这些专业参数。

原始代码如下：

import numpy as np
import matplotlib.pyplot as plt

# sampling rate
dt = 0.01  # 0.01초 간격

# rpm raw data
data1 = np.array([735, 743, 744, 739, 734, 752, 750, 745, 739, 751, 751, 745, 740, 757, 755, 749, 743, 745, 748, 742, 736, 748, 747, 741, 734, 751, 752, 747, 740, 737, 734, 730, 722, 731, 737, 732, 723, 738, 747, 743, 734, 747, 751, 749, 738, 750, 759, 755, 746, 752, 757, 753, 744, 747, 756, 753, 743, 749, 754, 752, 742, 743, 750, 747, 738, 740, 745, 745, 734, 737, 746, 746,734, 747, 754, 752, 741, 746, 753, 753, 741, 745, 752, 751, 741, 749, 759, 758, 746, 754,761, 760, 748, 753, 761, 761, 750, 751, 756, 757])


# FFT
fft_result1 = np.fft.fft(data1)
fft_magnitude1 = np.abs(fft_result1)

N = len(data1)
freq = np.fft.fftfreq(N, d=dt)  

# plot
plt.figure(figsize=(10, 6))
plt.plot(freq, fft_magnitude1, label='Signal 1')
plt.title('FFT Result')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude')
plt.xlim(4, 50)  
plt.ylim(0, 500)
plt.xticks(np.arange(4, 51, 2))  
plt.legend()
plt.grid(True)
plt.show()

解决方案 #

使用scipy.signal库实现所有商用级参数配置，同时完成幅度单位转换。以下是完整实现代码及关键步骤说明：

完整实现代码 #

import numpy as np
import matplotlib.pyplot as plt
from scipy import signal

# 原始采样参数
original_dt = 0.01
original_fs = 1 / original_dt  # 原始采样率100Hz
target_fs = 2048  # 目标重采样频率2048Hz

# RPM原始数据
data1 = np.array([735, 743, 744, 739, 734, 752, 750, 745, 739, 751, 751, 745, 740, 757, 755, 749, 743, 745, 748, 742, 736, 748, 747, 741, 734, 751, 752, 747, 740, 737, 734, 730, 722, 731, 737, 732, 723, 738, 747, 743, 734, 747, 751, 749, 738, 750, 759, 755, 746, 752, 757, 753, 744, 747, 756, 753, 743, 749, 754, 752, 742, 743, 750, 747, 738, 740, 745, 745, 734, 737, 746, 746,734, 747, 754, 752, 741, 746, 753, 753, 741, 745, 752, 751, 741, 749, 759, 758, 746, 754,761, 760, 748, 753, 761, 761, 750, 751, 756, 757])

# 1. 重采样到2048Hz
resampled_len = int(len(data1) * target_fs / original_fs)
resampled_data = signal.resample(data1, resampled_len)

# 2. 配置频谱分析参数
nfft = 8192  # 频谱大小
window = signal.windows.hann(nfft)  # 汉宁窗
noverlap = int(nfft * 0.5)  # 50%重叠

# 3. 计算频谱（带窗、重叠处理）
freqs, psd = signal.welch(
    resampled_data,
    fs=target_fs,
    window=window,
    nperseg=nfft,
    noverlap=noverlap,
    nfft=nfft,
    scaling='spectrum'  # 输出频谱幅度而非功率
)

# 4. 转换为峰值幅度并修正单位为[n/rpm]
window_sum = np.sum(window)
# 修正窗函数影响，转换为峰值幅度
peak_magnitude = (np.sqrt(psd) * window_sum) / target_fs * 2

# 5. 绘图
plt.figure(figsize=(12, 6))
plt.plot(freqs, peak_magnitude, label='RPM Signal Spectrum')
plt.title('RPM Signal FFT with Commercial Parameters')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude [n/rpm]')
plt.xlim(4, 50)
plt.xticks(np.arange(4, 51, 2))
plt.legend()
plt.grid(True)
plt.show()

关键参数实现说明 #

• 重采样2048Hz：使用signal.resample将原始100Hz信号插值到目标采样率，保证频谱分析的基准一致性。
• 汉宁窗：通过signal.windows.hann生成窗函数，抑制频谱泄漏，提升频率分辨率。
• 50%重叠：设置noverlap为频谱大小的一半，让相邻分析段覆盖50%的数据，提升频谱连续性。
• 频谱大小8192：指定nfft=8192，控制FFT计算长度，获得更高的频率精度。
• 峰值幅度缩放：通过窗函数能量修正，将Welch方法输出的频谱值转换为峰值幅度，匹配商用工具的显示逻辑。
• Y轴单位[n/rpm]：通过幅度修正公式，将FFT幅度值映射到RPM信号的物理单位，直观反映各频率下的信号强度。

内容的提问来源于stack exchange，提问作者wonsik lee