在 python 中实现 Cooley-Tukey 风格的 FFT 的问题

问题描述 投票:0回答:0

出于教育目的,我决定尝试在 python 中实现 Cooley-Tukey 风格的 FFT 算法,我必须说我绝不是专业的 python 开发人员,也不是信号处理专家,所以难看的代码是有点预料到的。

无论如何,我按照维基百科关于该主题的文章从原始DFT方程重新证明了自己在纸上的必要方程。

我尝试使用方程式和提供的伪代码在 Python 中实现它。然后我将一个正弦波注入我的 FFT 代码,而不是在香农频率的任何一侧获得我的正弦波频率,我得到这个频率加上另一个频率,它从我的信号中吸收更多的能量,我增加我的越多采样频率。

为了排除故障,我只将样本分成偶数和奇数索引一次。

这里是正弦波生成:

fs = 500
Ts = 1/fs
samples = 128
t = np.arange(0, samples * Ts, Ts)

f1 = 100
#f2 = 110

sine1 = np.array([np.sin(2 * np.pi * f1 * t[i]) for i in range(0,samples)])
#sine2 = np.array([np.sin(2 * np.pi * f2 * t[i]) for i in range(0,samples)])

signal = sine1

这是我的 FFT 实现:

#FFT algorithm

N = len(signal)

#pre-calculation check-ups
if N < 2:
    print("Selected sample is too small (less than 2 points)")
else:
    powerOfTwo = np.log2(N)
    if powerOfTwo.is_integer() == False:
        print("Selected sample isn't a power of 2")
    
    #start of the actual algorithm
    else:
        FFT1 = []
        FFT2 = []
        FFT = []
        evenArr = []
        evenArrFFT = []
        oddArr = []
        oddArrFFT = []
        TF = []
        
        #sorts samples by parity
        for i in np.arange(0, N):
            if(i % 2 == 0):
                evenArr = np.append(evenArr, signal[i])
            else:
                oddArr = np.append(oddArr, signal[i])
                
        #applies DFT to the even and odd parts of the original sample array individually
        evenArrFFT = DFT(evenArr, "blackman")
        oddArrFFT = DFT(oddArr, "blackman")
        
        #initializes the twiddle factor
        for k in np.arange(0, int(N/2)):
            TF = np.append(TF, np.exp(-1j * 2 * np.pi * k/N))
            
        FFT1 = evenArrFFT + TF * oddArrFFT  #second half of the FFT
        FFT2 = evenArrFFT - TF * oddArrFFT  #first half of the FFT
        FFT = np.append(FFT2, FFT1)         #complete FFT        

这是用我自己的 DFT 实现(带黑曼窗)得到的频谱:

DFT

这是用 Numpy 的 fft.fft(带黑曼窗口)获得的频谱:

然而,当我绘制我的 FFT 实现的输出时,这就是我得到的(再一次,你猜对了,用一个 blackman 窗口):

我确实尝试检查我的拆分为偶数和奇数索引是否有效,我的 DFT 实现是否有效,我是否在初始化旋转因子时没有出错,但我似乎无法找到问题的根源。

我完全知道互联网上存在更有效的 FFT 实现,本练习的目的是以我理解的方式自己编写代码。任何帮助将不胜感激。

python numpy signal-processing fft
© www.soinside.com 2019 - 2024. All rights reserved.