scipy.fft.fft函数与scipy.fft.rfft函数的异同

news2025/1/11 2:48:08

import numpy as np
from scipy import signal
import matplotlib.pyplot as plt
'''
思路：1）先利用fft计算得出其幅频值
     2）在利用rfft计算得出其幅频值，看1）和2）那个能还原出信号的原始幅值
'''
# 生成一个示例信号
np.random.seed(0)
fs = 512  # 采样频率
N = 512  # 信号长度
time = np.arange(N) / fs
sig = 10 * np.sin(2 * np.pi * 50 * time) + 5 * np.random.randn(N)

############ 1）利用fft计算得出其幅频值 ############
fft_signal = fft(sig)
fft_p = np.abs(fft_signal/N)
# 计算单边频谱
fft_p2 = fft_p[:N//2 + 1]
fft_p2[1:-1] = 2*fft_p2[1:-1]
fft_freqs = np.arange(0, N//2 + 1)*(fs / N)

############ 2）利用rfft计算得出其幅频值 ############
rfft_signal = rfft(sig)
rfft_p = np.abs(rfft_signal/N)
rfft_p2 = rfft_p
rfft_p2[1:-1] = 2*rfft_p2[1:-1]

plt.figure(1)
plt.subplot(211)
plt.plot(fft_freqs, fft_p2)
plt.title("fft method")
plt.subplot(212)
plt.plot(fft_freqs, rfft_p2)
plt.title("rfft method")

plt.show()

在这里插入图片描述
结论：利用fft方法和rfft方法皆可以反映出信号的真实幅值大小

############ 3）利用welch计算得出其功率谱 ############
f_welch, Pxx_welch = welch(sig, fs, nperseg=N)

############ 4）利用rfft计算得出其功率谱 ############
win = signal.get_window('hann', N)          # 获取窗函数
det_sig = detrend(sig, type='constant')     # 去除线性趋势
sig_win = win*det_sig                       # 加窗操作
rfft_win = rfft(sig_win)                    # 快速傅里叶变换，获取单边频谱

rfft_abs = np.conjugate(rfft_win)*rfft_win    # 模的平方(此处的rfft_abs为复数)
# rfft_abs = abs(rfft_win)**2

scale = 1.0 / (fs * (win*win).sum())          # 窗函数的缩放因子
Pxx_rfft = rfft_abs * scale                   # 乘于窗函数的缩放因子，保持加窗后的能量不变
Pxx_rfft[1:-1] = 2 * Pxx_rfft[1:-1]           # 获取双边频谱值
Pxx_rfft = Pxx_rfft.real                      # 获取实部

plt.figure(2)
plt.subplot(211)
plt.plot(f_welch, Pxx_welch)
plt.title('welch method')

plt.subplot(212)
plt.plot(f_welch, Pxx_rfft)
plt.title('rfft method')

plt.show()