傅里葉變換基礎(chǔ)
采樣定理：采樣頻率大于信號(hào)中最高頻率的2倍時(shí)白对，采樣后的數(shù)字信號(hào)可以完整地恢復(fù)出原始信號(hào)撩满。一般實(shí)際應(yīng)用中保證采樣頻率為信號(hào)最高頻率的2.56~4倍并鸵。
假設(shè)信號(hào)頻率為F担神，采樣頻率為Fs怕磨，采樣點(diǎn)數(shù)為N。為了方便進(jìn)行FFT運(yùn)算置尔，通常N取2的整數(shù)次冪。
N個(gè)采樣點(diǎn)氢伟，經(jīng)過(guò)fft變換后的結(jié)果為N個(gè)復(fù)數(shù)榜轿，每個(gè)復(fù)數(shù)對(duì)應(yīng)一個(gè)頻率(第n<=N/2個(gè)點(diǎn)對(duì)應(yīng)的頻率為(n-1)/N*Fs)，該復(fù)數(shù)的模值表示該頻率的振幅特征朵锣。該振幅特征和原始信號(hào)的振幅之間的關(guān)系是：如果原始信號(hào)的振幅為A谬盐，則fft結(jié)果的每個(gè)點(diǎn)(除了第一個(gè)直流分量點(diǎn))的模值就是A的N/2倍；而第一個(gè)點(diǎn)的模值是直流分量振幅的N倍诚些。
注意：這N個(gè)復(fù)數(shù)點(diǎn)去掉第一個(gè)點(diǎn)后剩下的N-1個(gè)點(diǎn)是關(guān)于其中心共軛對(duì)稱的飞傀，因此實(shí)際只需要取前一半點(diǎn)的頻譜即可，因?yàn)楣曹棇?duì)稱的兩個(gè)點(diǎn)的模值(振幅)相同诬烹。
In: np.fft.fft(np.array([1,2,3]))
Out: array([ 6. +0.j , -1.5+0.8660254j, -1.5-0.8660254j])
In: np.fft.fft(np.array([1,2,3,4]))
Out: array([10.+0.j, -2.+2.j, -2.+0.j, -2.-2.j])
numpy實(shí)現(xiàn)
我們先模擬一個(gè)一維時(shí)序信號(hào)y砸烦，它由2V的直流分量(0Hz)，和振幅為3V绞吁，頻率為50Hz的交流信號(hào)幢痘，以及振幅為1.5V，頻率為75Hz的交流信號(hào)組成：
y = 2 + 3*np.cos(2*np.pi*50*t) + 1.5*np.cos(2*np.pi*75*t)
然后我們采用256Hz的采樣頻率家破，總共采樣256個(gè)點(diǎn)颜说。
傅里葉變換結(jié)果如下圖所示，no normalization對(duì)應(yīng)的是原始的fft結(jié)果汰聋，normalization是將fft結(jié)果的振幅特征轉(zhuǎn)換為原始信號(hào)的振幅门粪，可以看到振幅為2V，3V烹困，1.5V的信號(hào)分別被解析了出來(lái)玄妈。
import numpy as np
import matplotlib.pyplot as plt
Fs = 256 # 采樣頻率 要大于信號(hào)頻率的兩倍可恢復(fù)信號(hào)
t = np.arange(0, 1, 1.0/Fs) # 1s采樣Fs個(gè)點(diǎn)
F1 = 50 # 信號(hào)1的頻率
F2 = 75 # 信號(hào)2的頻率
y = 2 + 3*np.cos(2*np.pi*F1*t) + 1.5*np.cos(2*np.pi*F2*t)
N = len(t) # 采樣點(diǎn)數(shù)
freq = np.arange(N) / N * Fs
Y1 = np.fft.fft(y) # 復(fù)數(shù)
print('Y[0]', Y1[0])
print('Y[128]', abs(Y1)[120:137])
Y = Y1 / (N/2) # 換算成實(shí)際的振幅
Y[0] = Y[0] / 2
freq_half = freq[range(int(N/2))]
Y_half = Y[range(int(N/2))]
fig, ax = plt.subplots(4, 1, figsize=(12, 12))
ax[0].plot(t, y)
ax[0].set_xlabel('Time (s)')
ax[0].set_ylabel('Amplitude')
ax[1].plot(freq, abs(Y1), 'r', label='no normalization')
ax[1].set_xlabel('Freq (Hz)')
ax[1].set_ylabel('Amplitude')
ax[1].legend()
ax[2].plot(freq, abs(Y), 'r', label='normalization')
ax[2].set_xlabel('Freq (Hz)')
ax[2].set_ylabel('Amplitude')
ax[2].set_yticks(np.arange(0, 3))
ax[2].legend()
ax[3].plot(freq_half, abs(Y_half), 'b', label='normalization')
ax[3].set_xlabel('Freq (Hz)')
ax[3].set_ylabel('Amplitude')
ax[3].set_yticks(np.arange(0, 3))
ax[3].legend()
# plt.show()
plt.savefig('a.png')
plt.close()
